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Abstract 



My thesis deals with a fundamental question of high energy gamma-ray astronomy. Namely, I 
studied the cut-off shape of the Crab pulsar spectrum to distinguish between the leading scenarios 
for the pulsar models. 

Pulsars are celestial objects, which emit periodic pulsed electromagnetic radiation (pulsa- 
tion) from radio to high energy gamma-rays. Two major scenarios evolved in past 40 years to 
explain the pulsation mechanism: the inner magnetosphere scenario and the outer magnetosphere 
scenario. Both scenarios predict a high energy cut-off in the gamma-ray energy spectrum, but 
with different cut-off sharpness. An exponential cut-off is expected for the outer magnetosphere 
scenario while a super-exponential cut-off is predicted for the inner magnetosphere scenario. 
Therefore, one of the best ways to confirm or rule out these scenarios is to measure the en- 
ergy spectrum of a pulsar at around the cut-off energy, i.e., at energies between a few GeV and 
a few tens of GeV. All past attempts (> 10) to measure pulsars with ground-based instruments 
have failed while satellite-borne detectors had a too small area to study detailed spectra in the 
GeV domain. 

In this thesis, the gamma-ray emission at around the cut-off energy from the Crab pulsar is 
studied with the MAGIC telescope. The public data of the satellite-borne gamma-ray detector, 
Fermi-LAI, are also analyzed in order to discuss the MAGIC observation results in comparison 
with the adjacent energy band. 

In late 2007, a new trigger system (SUM trigger system) allowed to reduce the threshold 
energy of the MAGIC telescope from 50 GeV to 25 GeV and the Crab pulsar was successfully 
detected during observations from October 2007 and January 2009. My analysis reveals that 
the energy spectrum is consistent with a simple power law between 25 GeV to 100 GeV. The 
extension of the energy spectrum up to 100 GeV rules out the inner magnetosphere scenario. 

Fermi-LAI started operation in August 2008. The Fermi-LAI data reveal that a power law 
with an exponential cut-off at a few GeV can well describe the energy spectrum of the Crab 
pulsar between 100 MeV and 30 GeV. This is consistent with the outer magnetosphere scenario 
and again, inconsistent with the inner magnetosphere scenario. 

The measurements of both experiments strongly disfavor the inner magnetosphere scenario. 
However, by combining the results of the two experiments, it turns out that even the standard 
outer magnetosphere scenario cannot explain the measurements. Various assumptions have been 
made to explain this discrepancy. By modifying the energy spectrum of the electrons which 
emit high energy gamma-rays via the curvature radiation, the combined measurements can be 
reproduced but further studies with higher statistics and a better energy resolution are needed to 
support this assumption. 

The energy-dependent pulse profile from 100 MeV to 100 GeV has also been studied in detail. 
Many interesting features have been found, among which the variabilities of the pulse edges and 
the pulse peak phases are the most remarkable. More data would allow an investigation of the 
fine structure of the pulsar magnetosphere based on these features. 

Aiming at better observations of pulsars and other sources below 100 GeV, a new photosen- 
sor, HPD R9792U-40, has been investigated. Many beneficial properties, such as a very high 
photodetection efficiency, an extremely low ion-feedback probability and an excellent charge 
resolution have been found. 
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Abstract 



Diese Dissertation behandelt eine grundlegende Fragestellung der Hochenergie-Astrophysik. 
Um herauszufinden, welches der fuhrenden Modelle die Emission von Pulsaren korrekt beschreibt, 
untersuche ich das Abknickverhalten des Energiespektrums des Pulsars im Krebsnebel zu hohen 
Energien hin. Pulsare sind astronomische Objekte, die durch periodische Abstrahlung ("Pulsa- 
tion") elektromagnetischer Strahlung von Radiowellen bis hin zu Gammastrahlen gekennzeich- 
net sind. Zwei Modellklassen haben sich in den vergangenen 40 Jahren herausgebildet, welche 
den Pulsationsmechanismus mit Teilchenbeschleunigung in der inneren bzw. aussssseren Mag- 
netosphare erklaren konnen. Beide Szenarien sagen ein Abknicken des Energiespektrums im 
Gammastrahlungsbereich voraus, allerdings mit unterschiedlicher Starke. Wenn die Strahlung in 
der ausssseren Magnetosphare erzeugt wird, erwartet man einen exponentiellen Abfall, wahrend 
fiir die Erzeugung in der inneren Magnetosphare ein noch starkerer Abfall vorhergesagt wird. 
Deswegen ist die Bestimmung eines Pulsar- Energiespektrums in der Nahe des erwarteten Abknick- 
ens, also bei Energien zwischen einigen GeV und einigen zehn GeV eine bevorzugte Methode, 
einem dieser Szenerien Glaubwiirdigkeit zu verleihen oder es andererseits auszuschliesssssen. 
Alle bisherigen Versuche, Pulsarspektren im genannten Energiebereich mit bodengebundenen 
Instrumenten zu vermessen, schlugen fehl; gleichzeitig hatten satellitengestiitzte Detektoren eine 
zu kleine Sammelflache, als dass die genaue Vermessung von Spektren moglich gewesen ware. 

Diese Arbeit untersucht die Gammastrahlungsemission des Krebsnebels in der Umgebung 
der Abknick-Energie mit Hilfe des MAGIC-Teleskops. Ebenfalls ausgewertet werden offentlich 
zugangliche Daten des satellitengestiitzten Gammastrahlungsdetektors Fermi-LAT, so dass die 
MAGIC-Ergebnisse mit dem bei niedrigeren Energien anschliessenden Energiebereich von Fermi- 
LAT verglichen werden konnen. 

Vor kurzem wurde ein neues Triggersystem (ein analoger Summentrigger) in Betrieb genom- 
men, welches die Energieschwelle des MAGIC-Teleskopes von 50 GeV auf 25 GeV herunter- 
setzt; damit wurde wahrend Beobachtungen von Oktober 2007 bis Januar 2009 der Pulsar im 
Krebsnebel erstmalig nachgewiesen. Meine Untersuchungen zeigen, dass das Energiespektum 
im Bereich von 25 GeV bis 100 GeV mit einem einfachen Potenzgesetz vertraglich ist. Allein 
schon die Tatsache, dass sich das Spektrum bis 100 GeV erstreckt, schliessssst eine Erzeugung 
in der inneren Magnetosphare aus. 

Die Fermi-LAT Beobachtungen begannen im August 2008. Die offentlichen Fermi-LAT 
Daten zeigen, dass ein Potenzgesetz mit exponentiellem Abfall bei einigen GeV das Krebspulsar- 
spektrum zwischen 100 MeV und 30 GeV gut beschreibt. Das ist vertraglich mit einer Erzeugung 
in der ausssseren, aber nicht in der inneren Magnetosphare. Die Messungen beider Instrumente 
deuten somit stark darauf hin, dass das entsprechende Modell fiir eine Erzeugung in der inneren 
Magnetosphare die Daten nicht korrekt beschreibt. Wenn man nun die Resultate beider Messun- 
gen kombiniert, zeigt sich, dass ein einfaches Modell zur Gammastrahlungs-Erzeugung in der 
ausssseren Magnetosphare die Daten nur unzulanglich beschreibt. Verschiedene Annahmen sind 
notwendig, um die vorhandenen Abweichungen zu erklaren. 

Weiterhin wurde die Energieabhangigkeit des PulspronTs zwischen 100 MeV und 100 GeV 
genau untersucht. Diverse interessante Eigenschaften wurden gefunden, von denen die Veranderungen 
in den Pulsflanken und die Entwicklung der Phasen der Pulsationsmaxima die erwahnenswertesten 
sind. Von diesen Eigenschaften ausgehend, wurde eine grossssere Menge von Beobachtungs- 
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daten erlauben, den Aufbau der Pulsarmagnetosphare genau zu untersuchen. 

Mit dem Ziel, die Beobachachtungsmoglichkeiten fur Pulsare und andere Quellklassen unter 
100 GeV zu verbessern, wurde eine neuer Photondetekor, der HPD R9792U-40, charakter- 
isiert. Viele positive Eigenschaften konnten gefunden werden, wie beispielsweise eine sehr hohe 
Photon-Nachweiseffizienz, eine ausssserordentlich niedrige Ionenriickkopplungswahrscheinlichkeit, 
sowie eine ausgezeichnete Ladungsaufiosung. 
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Introduction 



The first detection of the TeV gamma-rays from a celestial object, i.e. the discovery of the Crab 
Nebula by the Whipple telescope in 1989, opened a new field in astronomy, which is called very 
high energy (VHE) gamma-ray astronomy. In this thesis, VHE gamma-ray stands for photons 
with energies above 10 GeV. As of now approximately 100 VHE gamma-ray sources are known, 
which are categorized into several classes such as active galactic nuclei, supernova remnants, 
pulsar wind nebulae and gamma-ray binaries. Observations in VHE gamma-rays have given first 
insight into the nature of these extremely dynamical objects emitting non-thermal radiation. 

Pulsars are a class of celestial objects, which emit periodic pulsed radiation (pulsation) ex- 
tending from radio up to gamma-rays. VHE gamma-rays can serve as an important probe for the 
radiation mechanism of pulsars, too. Pulsars are explained as rapidly rotating neutron stars which 
possess extremely strong magnetic fields. Electrons are accelerated within their magnetosphere 
by strong electric field and emit beamed electromagnetic radiation, which will be observed as 
pulsation due to the rotation of the neutrons star. On top of this general picture of the pulsation 
mechanism, there are two competing major scenarios which specify the acceleration/emission 
region within the magnetosphere. One is the inner magnetosphere scenario, in which the pulsa- 
tion originates from near the magnetic pole on the neutron star surface. The other is the outer 
magnetosphere scenario, in which the pulsation comes from a region along the last closed mag- 
netic field lines in the outer magnetosphere. Both scenarios could reasonably explain all the 
features of pulsars observed before 2007. One of the best ways to verify or refute these scenarios 
is measuring the energy spectrum at around cut-off energy, i.e. at energies between a few GeV 
and a few tens of GeV. The reason is as follows: Both scenarios predict a high energy cut-off 
in the gamma-ray energy spectrum, but with different cut-off sharpness. An exponential cut-off 
is expected for the outer magnetosphere scenario while a super-exponential cut-off is predicted 
for the inner magnetosphere scenario. The sharpness of the cut-off and the highest energy of the 
observed photons allow one to constrain the emission region. 

Before 2007, a satellite-borne detector, EGRET, detected 7 gamma-ray pulsars and the en- 
ergy spectra could well be measured only up to ~ 5 GeV. On the other hand, ground-based 
imaging atmospheric Cherenkov telescopes (IACTs) could set flux upper limits only above 100 
GeV. There existed no sufficient measurement in the important energy range, i.e., at energies 
between a few GeV and a few tens of GeV. It was evident that filling this energy gap would lead 
to a clarification between the competing scenarios and a better understanding of the pulsation 
mechanism. 

The MAGIC telescope is the IACT that has a largest single dish reflector with a 17 m diam- 
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eter. Accordingly, it has the lowest energy threshold among IACTs (50 GeV in the case of the 
standard trigger). MAGIC has been the best instrument to fill the energy gap from the higher 
side. Moreover, in October 2007, the new trigger system (SUM trigger system) was imple- 
mented, which reduced the energy threshold even further, from 50 GeV to 25 GeV. Nearly at 
the same time, a new satellite-borne gamma-ray detector, Fermi-LAT, became operational in Au- 
gust 2008, which has a ~ 10 times better sensitivity than EGRET and could measure the pulsar 
energy spectrum well beyond 10 GeV. It began filling the gap from the lower side. 

In this thesis, the observational results on the Crab pulsar by MAGIC at energies above 25 
GeV are presented. The public data of the Fermi-LAT on the Crab pulsar are also analyzed 
from 100 MeV to ~ 30 GeV. The combined analysis of the results from two experiments is also 
carried out, carefully taking into account the systematic uncertainties of both experiments. Then, 
several constraints in the pulsar model based on the combined energy spectrum are discussed. In 
addition to the spectral study, energy-dependent pulse profiles between 100 MeV to 100 GeV are 
intensely studied. The possibility to infer the fine structure of the emission region based on the 
pulse profile is also discussed. 

Aiming for better observations of some selected pulsars and other sources below 100 GeV 
with MAGIC, a new photodetector, HPD R9792U-40, is investigated. HPD R9792U-40 is a hy- 
brid photodetector. Since a hybrid photodetector has never been used in any IACTs, its properties 
and performance are thoroughly studied. 

The thesis is structured as follows. An introduction to VHE gamma-ray astronomy is given 
in Chapter 1 . Pulsars and theoretical models of the pulsation mechanism are introduced in Chap- 
ter 2. Chapter 3 describes the IACT technique and the MAGIC telescope. The analysis methods 
of the MAGIC data and its performance are explained in detail in Chapter 4. The analysis results 
of the MAGIC data and the public Fermi-LAT data on the Crab pulsar are presented in Chapter 5 
and Chapter 6, respectively. In Chapter 7, a combined analysis of the results of the two experi- 
ments is performed. Physics discussions on the results are presented in Chapter 8. The properties 
and performance of HPD R9792U-40 and the development of associated operation circuits are 
presented in Chapter 9. Conclusions and outlook are added in Chapter 10. 



Chapter 1 

Very High Energy Gamma Ray Astronomy 

see http://wwwmagic.mppmu.mpg.de/publications/theses/TS aito .pdf 
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see http://wwwmagic.mppmu.mpg.de/publications/theses/TS aito .pdf 
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The MAGIC Telescope 
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Analysis Method of MAGIC Data 

see http://wwwmagic.mppmu.mpg.de/publications/theses/TS aito .pdf 
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Chapter 5 



MAGIC Observations of the Crab Pulsar 
and Data Analysis 

The MAGIC observations are grouped in one-year long cycles and Cycle I started in May 2005. 
In Cycles I and II, the Crab pulsar was observed with the standard trigger. 16 hours of good- 
condition data showed only a weak signal of pulsation with 2.9 a (see Sect ??). Here, I analyze 
the data recorded in Cycle III and Cycle IV with the SUM trigger with a much lower energy 
threshold and a higher sensitivity below 100 GeV. 

5.1 Observations 

The SUM trigger was installed in October 2007 and, subsequently, the Crab pulsar was observed 
in Cycles III and IV, for 48 hours (over 47 days) and 78 hours (over 36 days), respectively. 
All the observations were made in ON-mode, since the SUM trigger was designed for ON- 
mode observations. Although, as described in Sect ??, pulsar observations do not require OFF 
observations, in order to assure the validity of the analysis chain and the quality of the data sets 
by the Crab nebula emission, OFF observations were made in October and December 2007 for 
~ 10 hours. 

5.2 Sum Trigger Sub-patch Malfunction 

One and five out of the 36 Sum Trigger sub-patches were malfunctioning during the Crab pulsar 
observations in Cycle III and Cycle IV, respectively. The effect can be seen by plotting COG of 
images on the camera. The regions near the broken patches show a hole in COG distribution, as 
shown in the top left (Cycle III data) and the middle left (Cycle IV data) panel of Fig l5.ll When 
a variability of the Crab pulsar flux between the two cycles is discussed, this difference must 
be taken into account. By deactivating the broken patches in MC, the effect can be reproduced 
as shown in the top right (Cycle III MC) and the middle right (Cycle IV MC). The bottom left 
and the bottom right panel show the difference between Cycle III and Cycle IV in data and MC, 
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respectively. According to the MC, the differences in gamma-ray detection efficiency between 
the two cycles for SIZE 25-50, 50-100, 100-200 and 200-400 ph.e. are about 21%, 17%, 11% 
and 7%. This effect will be corrected when the variability is discussed (see Sect. I5.9I ). The 
calculation of the energy spectra will also be carried out with MCs which reflect these subpatch 
malfunctions. 



5.2 Sum Trigger Sub-patch Malfunction 
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Figure 5.1: COG distributions on the camera. The left and right columns show observed data and gamma- 
ray MC events, respectively. The first, second and third rows show Cycle 111, Cycle TV and difference 
between the two cycles, respectively. The effect of the subpatch malfunction is seen as holes of COG 
distributions. Simulations reproduce the effect reasonably well. 
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5.3 Data Selection 

The new SUM trigger system lowers the energy threshold from 50 GeV to 25 GeV (see Sect. 
??). It has a big impact on the pulsar observations because most of the signal is expected below 
50 GeV (see Fig. ??). However, such low energy events are easily affected by the observational 
conditions. A slightly worse condition may result in a significant increase in the energy threshold 
and a worse sensitivity. Therefore, a more careful data selection than the normal observations 
above 50 GeV is required, in terms of both hardware and environmental conditions. 



5.3.1 Reflector Performance Selection 

The mirror panels of the MAGIC telescope are adjustable (see Sect. ??) and their alignment on 
some days can be worse than it should be. As described in Sect. ??, this can be checked by muon 
events. Blue points in Fig. 15.21 show the light collection efficiency estimated by muon events 
(the conversion factor from the number of photons hitting the reflector to that of the detected 
photoelectrons) for each observation day. The days when the efficiency is lower than 0.081 (5% 
lower than the average) are excluded from the analysis, since it may affect the trigger efficiency 
and the energy reconstruction. On average, Cycle IV has 3% lower efficiency than Cycle III, 
which is taken into account in the MC when the energy spectrum is calculated. The PSF was 
checked by fitting a linear function to the ArcWidth/ Radius — Radius relation (see Sect ??). 
The value at Radius = 1.15 degrees is used to evaluate the PSF, as shown by the red points in 
Fig. 15.21 From the study with MC, it was found that the PSF in Cycle IV was ~ 5% worse than 
that in Cycle III, which is also taken into account when the energy spectrum is computed. 
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Figure 5.2: Daily monitoring of the light collection efficiency and the PSF of the reflector by muons. The 
horizontal axis indicates the observation day and the vertical axis indicates the mean ArcWidth and the 
mean conversion factor. Data taken in the days when the light collection efficiency (blue points) is below 
the limit (a blue line) are excluded from the analysis. 
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5.3.2 Zenith Angle Selection 

As the angular distance from the source to the zenith (Zenith Angle, ZA) becomes larger, the 
distance from the shower maximum to the telescope increases. The relation between the gamma- 
ray energy and the number of detected photons would change depending on the ZA, due to the 
different Cherenkov photon density on the ground Q. Needless to say, the threshold energy is 
also affected. Therefore, in order to assure the lowest possible threshold and a uniform SIZE - 
energy relation, I selected data with the ZA below 20 degrees. 
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Figure 5.3: ZA distribution of the data. Most of the data are taken below 20 degrees in ZA. In order to 
assure the lowest possible energy threshold and the uniform SIZE - energy relation, data with ZA below 
20 degrees are used in the analysis. 
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5.3.3 Cloudiness Selection 

As described in Chapter ??, haze and clouds affect the air shower images. Especially for the 
lowest energies below 50 GeV, the effect would be significant. Fig. 15 .4! shows the event rate H as 
a function of Cloudiness (see Sect. ??). A clear anti-correlation between Cloudiness and the 

'This is a consequence of three effects: a) a higher Cherenkov threshold leading to fewer photons (see the top 
right panel of Fig. ??), b) a wider spread of photons on the ground due to the larger distance and c) increase in 
absorption and scattering losses (see Sect. ??). 

2 Not the trigger rate but the rate of the events which are not completely erased by the image cleaning procedure. 
Therefore, most of the NSB accidental trigger events are not included in the rate 
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event rate can be seen. To assure that no data is affected by haze and clouds, I selected data taken 
with Cloudiness lower than 20, as shown in Fig. 15.51 If the daily average value is more than 20, 
all the data taken on that day are excluded. 



Event Rate vs Cloudiness 
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Figure 5.4: The event rate as a function of Cloudiness. A clear anti-correlation can be seen. Cut values 
in the event rate and in Cloudiness are indicated by the red dotted lines. 




Figure 5.5: Left: Cloudiness for each run. The horizontal axis indicates the observation day and the 
vertical axis indicates Cloudiness. A red dot corresponds to one run f~ 1 minute). The mean and the 
RMS of the Cloudiness for each day are denoted by blue crosses. The data with Cloudiness larger than 
20 are excluded from the analysis as denoted by the blue dotted line. Right: Cloudiness distribution. 
Most of the runs have cloudiness below 20 and the data with Cloudiness larger than 20 are excluded 
from the analysis. 



5.3 Data Selection 



17 



5.3.4 Event Rate Selection 

Even if the environmental conditions and the reflector status are good, the telescope performance 
might not be optimal because of inappropriate DAQ settings, such as wrong threshold settings 
and the change in signal transmission time (length), which can happen during long-term obser- 
vations. Such problems can be identified by checking the event rate. For the very low threshold 
observations like pulsar observations, almost half the trigger rate is due to NSB + after-pulsing 
accidental events and, hence, the trigger rate may not reflect improper DAQ settings. The event 
rate after the image cleaning would be more indicative because images of such accidental events 
would be completely erased by the image cleaning. As you can see from Fig. 15 .61 after cutting 
away bad reflector days and cloudy days, almost all runs have a good rate, meaning that hardware 
settings had been fine. Runs with the event rate lower than 380 Hz were discarded just in case. 



■n 700 p 
I 




Day Event Rate [Hz] 



Figure 5.6: Left: The event rate for each run. The horizontal axis indicates the observation day and the 
vertical axis indicates the event rate. A red dot corresponds to one run f~ 1 minute). The mean and the 
RMS of the event rate for each day are denoted by blue crosses. The data with the event rate smaller 
than 380 Hz are excluded from the analysis, as denoted by the blue dotted line. Right: The event rate 
distribution. Only a few runs have an event rate smaller than 380 Hz, and they are excluded from the 
analysis. 



5.3.5 Nebula Measurement Selection 

From the rate cut described in the previous section, it is almost guaranteed that event selection has 
been properly carried out. Just as an additional cross-check, the detection efficiency of the Crab 
nebula emission was examined. I analyzed the Crab nebula emissio n| with SIZE above 300, 
for which Hadronness cut is very powerful, i.e., the gamma-ray/hadron separation is highly 
efficient (Sect. ??) and analysis is rather easy. The Hadronness cut at 0.1 is applied to the data 
and the excess is evaluated with the ALPHA cut at 10 degrees. The background was estimated 
by fitting the ALPHA distribution from 20 degrees to 80 degrees with a parabolic function and 

3 As noted in Sect. ??, the nebula and the pulsar cannot be spatially resolved by IACTs. The non-pulsed gamma- 
ray emission above 100 GeV is considered as a nebula emission 
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5. MAGIC Observations of the Crab Pulsar and Data Analysis 



then extrapolating the fitted function down to degree. An example of this analysis is shown 
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Figure 5.7: An example of the signal from the Crab nebula with SIZE > 300 and Hadronness < 0.1. 
ALPHA analysis is used (see Sect. 11). 

in Fig. 15.71 which is for the observation on 1st January, 2008. In order to check the stability of 
detection efficiency, the significance of the excess for each of the observation days was scaled to 
one hour's observation and plotted in Fig. 15.81 They are very stable and their mean value is 16.1, 
which is consistent with the telescope sensitivity. 



I Crab Nebula Excess > 300 PhE | 




Figure 5.8: Left: Daily monitoring of the detection efficiency of the nebula emission. The horizontal axis 
indicates the observation day and the vertical axis indicates the excess significance normalized to one 
hour of observation time. Right: Distribution of the excess significance normalized to one hour. The 
variation is compatible with the statistical fluctuations. 

5.3.6 Summary of the Data Selection 

The selection of the good observation days are summarized in Table. 15.11 26 days out of 47 days 
and 18 days out of 36 days passed all the selection cuts for Cycle III and IV, respectively. Within 
a good day, some runs were also excluded due to unfavorable ZA, Cloudiness or the event rate. 
25.1 and 34.0 hours of data remained for Cycle III and IV, respectively. 



5.3 Data Selection 
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Date 


PSF/ 










Date 


PSF/ 








yymmdd 


Ref 


Z.A. 


Cloud 


Used 




yymmdd 


Ref. 


Z.A. 


Cloud 


Used 


Cycle ID 






080129 








yes 


071021 








yes 




080131 






Cloudy 


no 


071022 








yes 




080201 






Cloudy 


no 


071105 






Cloudy 


no 




080202 








yes 


071106 






Cloudy 


no 




080203 








yes 


071107 






Cloudy 


no 




Cycle IV 


071108 


Bad 






no 




081105 


Bad 






no 


071109 






Cloudy 


no 




081107 


Bad 






no 


071110 


Bad 






no 




081126 


Bad 






no 


071111 






Cloudy 


no 




081129 


Bad 






no 


071112 






Cloudy 


no 




081130 


Bad 






no 


071113 






Cloudy 


no 




081201 


Bad 






no 


071114 






Cloudy 


no 




081202 


Bad 






no 


071115 








yes 




081203 


Bad 






no 


071116 








yes 




081204 








yes 


071117 






Cloudy 


no 




081205 








yes 


071118 






Cloudy 


no 




081206 






Cloudy 


no 


071205 








yes 




081207 






Cloudy 


no 


071206 








yes 




081208 




> 20° 




no 


071207 








yes 




081210 




> 20° 




no 


071208 








yes 




081219 








yes 


071209 








yes 




081220 








yes 


071210 






Cloudy 


no 




081222 








yes 


071211 






Cloudy 


no 




081223 






Cloudy 


no 


071212 






Cloudy 


no 




081229 








yes 


071213 


Bad 






no 




090101 






Cloudy 


no 


071214 








yes 




090102 








yes 


071230 








yes 




090103 








yes 


071231 








yes 




090104 








yes 


080101 








yes 




090118 


Bad 






no 


080102 








yes 




090119 








yes 


080103 








yes 




090120 








yes 


080104 








yes 




090121 








yes 


080105 








yes 




090122 








yes 


080106 








yes 




090124 








yes 


080107 






Cloudy 


no 




090125 








yes 


080108 








yes 




090126 








yes 


080109 






Cloudy 


no 




090127 






Cloudy 


no 


080110 








yes 




090128 








yes 


080111 








yes 




090130 








yes 


080112 








yes 




090131 




> 20° 




no 


080113 








yes 




090201 




> 20° 




no 


080126 






Cloudy 


no 















Table 5.1: Selection of the data on daily basis based on the reflector performance, the zenith angle and 
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5. MAGIC Observations of the Crab Pulsar and Data Analysis 



5.4 ALPHA Cut Optimization for Pulsar Analysis 

For the pulsar analysis, the signal extraction can be carried out using the light curve (see Sect. 
??). The image parameter ALPHA is additionally used as an event selection parameter in order 
to increase the signal-to-noise ratio. The best cut value on ALPHA changes with SIZE because 
the larger the SIZE, the better the shower direction estimation. The dependency is especially 
strong at SIZE < 100, as one can see from the top right panel of Fig. 15.91 Since most of 
the signal from the Crab pulsar is expected at SIZE < 100, a S'/Z.E'-dependent ALPHA cut is 
applied to the Crab pulsar data. It is optimized as follows: First, MC gamma-ray and data (mostly 
hadron) samples are divided into 20 of log w (S I Z E) bins from 1 (10 ph.e.) to 3 (1000 ph.e.), as 
shown in Fig. 15.91 For each of the bins, the best ALPHA cut is calculated which maximizes 
the so-called Q-factor Q = e g /^/eh, where e g and e h is the fraction of events which survive 
the ALPHA cut for gamma-ray and hadron samples, respectively. Red stars in the top right 
panel of Fig. 15 .91 indicate the best ALPHA cut values for each bin. Then, those best values as a 
function of \og w (SIZE) are fitted by a function A cut (SIZE) = a(log 10 (5 1 Z E) +b) c , obtaining 
a = 3.7 x 10 4 , b = 1.674, c = —5.988 as the best parameters. The function is shown in the same 
panel as a black line. In the data analysis, the events which fulfill ALPHA < A cut (S I Z E) are 
considered as gamma-ray candidates. At size 25 ph.e. the cut is at 45 degrees, while at 250 ph.e. 
it is at 8 degrees. e h , e 7 and Q are shown in the bottom panel of Fig. 15.91 Q is approximately 1.5 
at 100 ph.e. and lower for smaller SIZEs. 

It should be noted that a cut in Hadronness was found not to improve the signal-to-noise 
ratio as much as in ALPHA at SIZE < 100 where most of the signal from the Crab pulsar is 
expected. In order to avoid systematic errors in analysis, a Hadronness cut is not applied. 



5.4 ALPHA Cut Optimization for Pulsar Analysis 
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Figure 5.9: Top left: ALPHA as a function of logio(SIZE) for an observed data sample. Since the 
data sample is dominated by the hadron background events, ALPHA distributes uniformly from to 
90 degree. A black line indicates the SI Z E-dependent cut values. Top right: ALPHA as a function 
of logio(SIZE) for an MC gamma-ray event sample. ALPHA is concentrated around 0, whereas the 
concentration gets weaker as SIZE decreases. Red stars indicate the best cut values which maximize the 
Q-factor. A black line indicates the SI Z E-dependent cut values. Bottom: The Q-factor (blue) and the 
survival ratios for hadrons (green) and gamma-rays (red) as a function of logio(SIZE). Q-factor is about 
1.5 at SIZE = 100 and 1 at SIZE = 25. 
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5. MAGIC Observations of the Crab Pulsar and Data Analysis 



5.5 Analysis of the Energy Spectrum of the Crab Nebula 

The Crab nebula is generally used as a standard candle for the calibration of the IACTs (see 
Sect. ??). It is the brightest steady point-like source above 100 GeV. Actually that is why it 
could conveniently be used to verify the quality of the data for the Crab pulsar analysis (see Sect. 
15.3 .5b . It would also be very important to calculate its energy spectrum in order to assure that the 
analysis tool and the MC samples are appropriate before starting the Crab pulsar analysis. 

5.5.1 Data Sample for the Crab Nebula Analysis 

I used 17.1 hours of data (ON-data) and 5.2 hours of OFF observation data (OFF-data), both of 
which were taken in October and December 2007 and passed all the selection criteria described 
in Sect. 15 .31 The conditions of the observations are summarized in Table 15 .21 The pointing 
position of the OFF-data is on the same declination as for the ON-data but differs by 2 hours 
in right ascension, which results in the same trajectory on the sky between ON and OFF with a 
2-hour time lag. 





ON 


OFF 


Dec. [deg.] 


22.014 


22.014 


R.A. [hour] 


5.5756 


3.5756 


Period 


October 2007 


October 2007 




and December 2007 


and December 2007 


Eff. Time [hours] 


17.1 


5.2 


Zenith Range [deg] 


6 to 20 


6 to 20 



Table 5.2: The observation conditions for the Crab nebula data. 



5.5.2 Energy Range of the Crab Nebula Analysis 

For the energy range below 50 GeV, a precise background estimation is not easy in the case 
of the nebula observations. The reason is as follows: For the nebula analysis, the light curve 
cannot be used for the signal extraction and, hence, the ALPHA distribution is used instead 
(see Sect. ??). However, since the images with a small SIZE are easily affected by the sky 
conditions, ALPHA distribution may slightly vary depending on the conditions, which may 
lead to a difference between ON-data and OFF-data. Since the observed events are strongly 
dominated by cosmic-ray background events, even a tiny mismatch in ALPHA distributions 
may result in a large systematic error in the background estimation. It should be stressed that 
this is not the case for the pulsar analysis. Since the light curve can be used for the pulsed 
signal extraction, it does not require any OFF-observation. Unless the sky condition changes 
significantly in the time scale of a pulsar period (34 ms for the Crab pulsar), which is very hard 
to imagine, the background estimation can be properly carried out. For this reason, the nebula 
emission was analyzed only above 50 GeV. 



5.5 Analysis of the Energy Spectrum of the Crab Nebula 
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The MC samples were generated from 6 GeV to 2 TeV. They are meant to be used for the 
analysis of the Crab pulsar, which is known not to have significant emission above 100 GeV. 
With these MC samples, it is very hard to analyze the spectrum above 1 TeV properly because no 
information above 2 TeV is present in MC and energy resolution is limited (~ 20% at 1 TeV). In 
order to avoid any analytical bias due to the limited energy range in MC production, the spectrum 
was analyzed up to 700 GeV 

The results are shown in Fig. 15.101 The spectrum was unfolded using the Tikhonov regular- 
ization method. It is known that the spectrum in this energy range is better fitted by a variable 
power law 

/ (E/300GeV)[ ri+r2lo ^^/ 300GeV )] (5.1) 



dEdAdt 
than a simple power law. 

d 3 F 



dEdAdt 



/o(£/300GeV) r (5.2) 



The power law fitting gives the best parameters of fo = (6.6 ± 0.7) x 10~ 10 [cm~ 2 s _1 TeV~ 1 ] 
and T = —2.11 ± 0.10, while x 2 = 4.74 with the degree of freedom 4 (x 2 probability = 31.5 
%). On the other hand, the variable power law gives f = (7.0 ± 0.8) x 10~ 10 [cm~ 2 s~ 1 TeV _1 ], 
r\ = -2.21 ± 0.15 and T 2 = -0.45 ± 0.47, while \ 2 = 3.55 with the degree of freedom 3 (x 2 
probability = 31.4%). In this limited energy range, both functions fit well. It should be noted 
that the fits take into account the correlation between the spectral points which is introduced by 
the unfolding procedure. The spectrum measured by HESS telescopes above 500 GeV (see 0T]) 
and the previous MAGIC measurements above 60 GeV (see ||2~0~I ) are also shown in the same 
figure. They are consistent with one another in the overlapping range, verifying the validity of 
the MC samples and analysis tools used for the pulsar analysis. 
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Figure 5.10: The energy spectrum of the Crab nebula. Red triangles, open squares and filled squares 
indicate the MAGIC observation results with the SUM trigger (this work), the published MAGIC result 
with the standard trigger taken in 2006 (taken from H20\l ), and the published HESS results (taken from 
HI1\I ). They are consistent in the overlapped energy region. The power law fit (a green line) and the 
variable power law fit (a blue line) to the MAGIC results with the SUM trigger are also shown. 



5.6 Analysis of Optical Pulsation 
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5.6 Analysis of Optical Pulsation 



For the pulsar analysis, the proper barycentric correction and the accurate pulsar period informa- 
tion are essential. In order to check if the barycentric correction and the pulse phase calculation 
are correctly done, the optical pulsation of the Crab pulsar is very useful and helpful since it can 
be clearly detected within 10 minutes of observation thanks to the large reflector of MAGIC. 



Optical Pulsation 



3 

2j 

£ 608 
c 

3 
O 

O 606 

Q 

< 

IL. 

604 
602 
600 
598 



if;;;; 




: ; ,1 ! 

; 1 ; \ 




» \ 








_ i 1 if ; ; 
- . |, i H f . j i i # 


J 


mm! "WWrtW^W 




— .... i .... i .. . , i . , , , i , , , . i , , . . 


, , , , i , , , , f , , , , f , , , , 



-0.2 



-0.1 



0.1 



0.2 



0.3 



0.4 



0.5 



0.6 



0.7 0.8 
Phase 



Figure 5.11: The light curve of the optical pulsation from the Crab pulsar observed by MAGIC with the 
central pixel (see Sect. 11). ~ 30 hours of the observation data are used. 10 minutes is enough to see the 
two peaks clearly. 



Fig. 15.1 II shows optical pulsation detected by the central pixel of MAGIC (see Sect. ??). 
~ 30 hours of the observation data from both cycles are used. The barycentric correction was 
done by Tempo and the period information i.e. v, v and to, are taken from the "Jodrell Bank 
Crab Pulsar Monthly Ephemeris", as described in Sect. ??. The phase is calculated by Eq. 
??. A glitch (see ??) occurred in May 2008, which is between the Cycle III and the Cycle IV 
observations. Since there is no Crab pulsar observation in that month and v, v and t are updated 
monthly, the glitch does not affect the phase calculations. 

The peak phase of PI is slightly shifted earlier with respect to the radio peak phase by ~ 0.01 
in phase, corresponding to ~ 300 fis, which is known and consistent with other observations 
(see e.g. H14610 . Pulse shapes are also in good agreement with other observations, although a 
small time variability has been reported and quantitative comparison is not easy (see e.g. HI 1011 ). 
Exactly the same method of phase calculation is applied to gamma-ray signals. 
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5.7 Detection of the Very High Energy Gamma-ray Pulsation 
from the Crab Pulsar 

After all the selection described in Sect. 15.31 59.1 hours (Cycle III + IV) of good data remained. 
The quality of the data sets, MC samples, analysis tools and the pulsar phase calculations were 
verified, as described in the previous sections. In this section, the analysis of the pulsed gamma- 
ray signal from the Crab pulsar is described. 



5.7.1 The Pulsed Gamma-ray Signal 




Phasogram : Cycle All (25 < Size < 500) 

x10 J ~ 



746 




: OP P1 

P1:6167 +-1444 (S 4.28 

P2: 11315 +- 1532 (S 7.: 

P12: 17482 +-2329 (S 7 




Figure 5.12: The light curves of the Crab pulsar detected by MAGIC. 59.1 hours of the observation data 
are used. Following the convention in pulsar analyses, the light curves are plotted from phase -1 to 1. 
The upper and lower panels show the light curves with 11 bins (3 milliseconds per bin) and with 22 bins 
(1.5 milliseconds per bin) per period. The red, shaded regions indicate the phases of PI and P2, the blue 
shaded regions indicate the Bridge emission phases and the black shaded regions indicates the OP ( off- 
pulse) phases (see Sect. 11 for the definition of these phases). The statistical significance of the excess for 
PI, P2, and PI + P2 are 4.28 a, 7.39 a and 7.51 a, respectively. 



5.7 Detection of the Very High Energy Gamma-ray Pulsation from the Crab Pulsar 
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The light curve of all events with SIZE range from 25 to 500 are shown in Fig. 15.121 A 
S I Z E-dependent ALPHA cut described in 15.41 was also applied. The background was esti- 
mated using the OP (off-pulse) region (phase 0.52 to 0.88). PI (phase -0.06 to 0.04), P2 (phase 
0.32 to 0.43) and the sum of PI and P2 have 6267±1444, 1 1315 ± 1532 and 17482 ±2329 excess 
events with statistical significance of 4.28cr , 7.39a and 7.51a, respectively. For the definition of 
the phase names, see Sect. ??. The flux of P2 is twice as high as that of PI. As can be seen in 
Fig. ??, at 1 GeV, PI has a higher flux than P2. The energy dependence of P2/P1 ratio will be 
discussed in Sect. 17.21 The so-called bridge emission, which is seen in some energy bands, is not 
visible in the MAGIC data. This will also be discussed in Sect. 17.21 As one can see the bottom 
panel of Fig. I5.12[ although PI is conventionally defined as 0.32 to 0.43, most of the excess is 
concentrated in a narrower phase interval. The precise discussion of the pulse shape will take 
place in Sect. 17.31 



5.7.2 Further Investigation of the Signal 

Since this is the first detection of a gamma-ray pulsar by an IACT0, it is important to assure that 
the signal is not an artifact of the analysis or of the instrument. A useful check is the growth 
of the number of excess events as a function of the number of background events. Since pulsar 
emission is thought to be stable in time, the excess should grow linearly. The results are shown 
in the left panel of Fig. 15.131 and indeed the excess grows linearly. The growth of statistical 
significance is also shown in the right panel. 




Figure 5.13: Left: the growth of the number of excess events as a function of the number of background 
events. Red, blue and black points indicate PI, P2, and PI + P2, respectively. Right: The growth of 
the statistical significance as a function of the number of background events. Red, blue, and black points 
indicate PI, P2, and PI + P2, respectively. 

4 The discovery of the Crab pulsar with MAGIC was first achieved by the efforts of my colleagues listed in Sect. 
?? and reported in J22) 
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5. MAGIC Observations of the Crab Pulsar and Data Analysis 



Another useful check is the "inverse selection" of events. If the excess is due to gamma-ray 
signals, events discarded by the ALPHA cuts should not contain a significant excess. Figure 
15.141 shows the light curve produced with events which do NOT pass the ALPHA cut. Excesses 
are compatible with the background fluctuation. 



Phasogram : Cycle All ( 25 < SIZE < 500 ) Anti-Alpha cut 

xio 3 




Figure 5.14: The light curve of the Crab pulsar for data which do not pass the ALPHA cut. No pulsed 
signal is visible, as expected. 



5.8 SIZE Dependence of the Pulsation 

In order to roughly estimate from which energies the excess events come from, the energy depen- 
dence of the light curves should be examined. Instead of the energy reconstructed with Random 
Forest (see Sect. ??), which suffers a trigger-bias effect at these low energies as discussed in 
Sect. ??, I use SIZE for this study. SIZE is the total number of photoelectrons in an image 
and a good indicator of the primary gamma-ray energy, especially for this data set, for which the 
ZA range is limited up to 20 degrees. As one can see in Fig. ??, SIZE in ph.e. corresponds 
roughly to two times the energy in GeV. 

The data with SIZE from 25 to 800 were divided into 5 bins in \og 10 (SIZE). Data with 
SIZE above 800 were also analyzed. The light curves of these six sub-samples are shown in 
Fig. 15.151 The numbers of excess events for PI, P2 and PI + P2, shown in the right upper corner 
of each panel, were calculated by estimating the background level with the OP region (0.52 to 
0.58). Most of the excess events are in the two lowest SIZE bins. The third (100 -200) and 
fourth (200-400) bins also show a 2 a level excess. Above 400 ph.e., no more excess is visible. 
The size-dependence of the excess is shown graphically in Fig 15 .161 

A detail calculation of the energy spectrum of the pulsed gamma-ray signal will take place in 
Sect. \5M 



5.8 SIZE Dependence of the Pulsation 
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Figure 5.16: The numbers of excess events as a function of SIZE. Red, blue and black lines indicate PI, 
P2 and PI + P2, respectively. Most of the excess events are at SIZE < 100, while some excess events 
are also visible at SIZE > 100. 



5.9 Time Variation of the Pulsation from the Crab Pulsar 
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5.9 Time Variation of the Pulsation from the Crab Pulsar 

The previous study by L. Kuiper et al. (see HI 1511 s ) shows the flux of the Crab pulsar at 1-10 MeV 
is stable in a time scale of years. However, it is still interesting to check if the flux measured 
by MAGIC is stable or not, especially because it is the flux beyond the cut-off energy. The 
growth of the number of excess events as a function of the number of background events already 
showed that there is no drastic change in flux (see Fig. I5.131 ). Here I examine the stability of 
the flux and the light curve more quantitatively. Due to the limited statistical significance of the 
excess, dividing data sample into too many subsets does not make much sense. Therefore, I only 
compared Cycle III and Cycle IV to search for a possible yearly variability. 



5.9.1 Variability in Light Curve 

Fig. 15. 171 shows the light curves of Cycle III (top) and Cycle IV (bottom). The SIZE range is 
from 25 to 500. In order to evaluate the variation in the light curve, the %2 test was performed for 
the 1 1 bins starting from -0.0682 to 0.432 in the histograms, which is roughly from the beginning 
of PI to the end of P2. The %2 was 5.00 while the degree of freedom was 10, indicating no 
significant difference between the two light curves. 



5.9.2 Variability in Flux 

As one can see from Fig. 15.171 although Cycle IV has 30% longer observation time, it shows 
less excess events than Cycle III. This can be explained only by statistical fluctuation but the 
hardware malfunction described in Sect. I5.2l may also have played a role. S'/ZE'-dependence of 
the number of excess events are also compared in Fig. 15.181 Left, middle and right panels are for 
PI, P2 and PI + P2. The difference in observation time between the two cycles is corrected. The 
effect of broken sub-patches estimated by MC (see Sect. 15.21 ) is also corrected, such that Cycle 

III and IV have the same gamma-ray detection efficiency. %2s are 1.04, 3.14, 2.46 for PI, P2 and 
PI + P2, while the number of dof is, 4 for all phase intervals. Signals observed in Cycle III and 

IV are statistically consistent and no significant variability is seen between Cycle III and IV. 



5.10 Energy Spectra of the Pulsation from the Crab Pulsar 

Here I show the energy spectra for the pulsation of PI, P2 and PI + P2. Since the effective area 
is increasing rapidly from 20 GeV to 200 GeV (see Sect. ??) and the energy resolution is rather 
poor in this energy region (see Sect. ??), the spectra must be calculated with great care. Events 
with SIZE lower than 30 will be excluded from the spectrum calculation in order to avoid a 
possible mismatch between MC and data (see Sect. ??). A few different methods will be tried in 
order to estimate the analytical uncertainty of the results. 
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Phasogram : Cycle III (25 < Size < 500) 
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Figure 5. 1 7: The comparison of the light curves between Cycle III and TV. Significant variation of the light 
curve cannot be seen. 



5.10.1 Excess Distribution in the Reconstructed Energy and Size 

Fig. 15 .191 shows the distributions of the number of excess events as a function of SIZE and the 
reconstructed energy. It should be emphasized again that the reconstructed energy is strongly 
biased by the trigger effect (see Sect. ??). The energy spectrum will be calculated from these 
distributions. SIZE (total charge in a shower image) has less systematic uncertainty than the 
reconstructed energy, whereas the resolution of the energy estimation by SIZE is poorer than 
that by the reconstructed energy. Therefore, it is important to analyze the spectrum with both and 
compare the results. 



5.10.2 Spectrum Calculation by Forward Unfolding Assuming a Power 
Law 

As long as the assumption of the spectral shape is valid, the forward unfolding method should 
provide the most robust and reliable result, as described in Sect. ??. Here, I assumed a power 
law spectrum with two free parameters, namely, differential flux at 30 GeV / 30 and the spectral 
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Figure 5.18: The comparison of the number of excess events as a function of SIZE between Cycle III and 
IV. The correction for the effect of the subpatch malfunction described in 15.21 is applied for each SIZE 
bin. In addition, the numbers of excess events for Cycle IV are multiplied by 25.1/34.0 in order to correct 
for the difference of the observation time. 



index T: 

^^UE/WGeVr (5.3) 

Due to the relatively poor statistics and the narrow energy range, even if the true spectrum is 
slightly curved, the power law assumption should be valid. 

The results based on the SIZE distribution (the left panel of Fig. 15.191) and the reconstructed 
energy distribution (the right panel of Fig. I5.191 l are shown in Table 1531 The obtained spectra 
are also graphically shown in Fig. I5.20l by green (based on SIZE) and red (based on the re- 
constructed energy) lines. Unfolded excess distributions, i.e., the expected excess distributions 
from the obtained spectra, are overlaid in Fig. 15.191 The \ 2 values between the unfolded excess 
distribution and the observed excess distribution are shown in the fifth column of Table [531 They 
are sufficiently small, assuring the validity of the power law assumption. 



5.10.3 Spectrum Calculation by Unfolding with Different Regularization 
Methods 

(Backward) unfolding does not require an assumption of the spectral shape a priori. However, 
as described in Sect ??, different regularization methods might give different results. The results 
are reliable only when all the regularization methods show consistency. Tikhonov (see H182H ). 
Schmelling (see H163II ) and Bertero (see [|42l) regularization methods were used for unfolding 
the reconstructed energy distribution. The results are shown in Fig 15.201 as black (Tikhonov), 
blue (Schmelling) and pink (Bertero) points. All the results are consistent. The results with the 
Tikhonov method were fitted by a power law function (Eq. I5.3I ). shown as a black line in the 
figure. It should be noted that points are correlated with each other because of the unfolding 
procedure but that the correlation is taken into account when fitting is performed. The best fit 
parameters and x 2 values are summarized in Table 15.31 
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Figure 5.19: The numbers of excess events as a function of SIZE (the left panel) or the reconstructed 
energy (the right panel). Points indicate the observed data used for the spectrum calculations, while solid 
lines indicate the expected number of excess events from the spectra obtained by the forward unfolding. In 
order to reduce the possible systematic error for the spectrum, the events with SIZE < 30 are excluded. 
Red, blue and black lines are for PI, P2 and PI + P2, respectively. 

5.10.4 Discussion of the Results 

As summarized in Table [5731 the three methods, namely, the forward unfolding with the SIZE 
distribution, the forward unfolding with the reconstructed energy distribution and the unfolding 
with reconstructed energy distribution using the Tikhonov regularization method, show consis- 
tent results. The energy spectra can be well described by a simple power law, which is partially 
due to the limited statistics and the narrow energy range. The flux of P2 at 30 GeV is twice as 
large as that of PI , being consistent with the light curve. On the other hand, no difference in spec- 
tral indices between them is visible. The indices are — 3.35±0.52, — 3.07±1.04 and — 3.45±0.54 
for PI + P2, PI and P2, respectively, in the case of the unfolding with the Tikhonov regularization 
method. 



5.11 Concluding Remarks 

After many careful checks, 59.1 hours of high quality data have been obtained. All the analy- 
sis tools have been carefully examined as well. From them, 6267 ± 1444, 11315 ± 1532 and 
17482 ± 2329 gamma-ray signal events have been detected for PI, P2 and PI + P2, respectively, 
corresponding to 4.28<r, 7.39<r and 7.51er in statistical significance. The light curves show some 
interesting features compared to lower energy bands, such as very narrow PI peak and the ab- 
sence of the bridge emission. These features will be further discussed in Chapter [7] together with 
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Phase 


Method 


/ 30 [lO^cnrt^TeV- 1 ] 


r 


X 2 /dof (prob.) 


PI +P2 


Forward with Size 
Forward with E rec 
Tikhonov 


17.3 ±2.1 

18.8 ±2.4 

14.9 ± 2.9 


-3.53 ± 0.41 
-3.42 ± 0.34 
-3.35 ± 0.52 


3.7/2 (15.7%) 
5.1/4(27.7%) 
2.90/1 (8.9%) 


PI 


Forward with Size 
Forward with E. rec 
Tikhonov 


5.7 ± 1.4 
6.0 ± 1.5 
4.5 ± 2.3 


-3.67 ± 0.80 
-3.06 ± 0.59 
-3.07 ± 1.04 


2.1/2(35.0%) 
7.2/4(12.6%) 
1.20/1 (27.3%) 


P2 


Forward with Size 
Forward with E rec 
Tikhonov 


11.3 ± 1.5 
12.6 ± 1.6 
10.0 ± 1.9 


-3.53 ± 0.39 
-3.54 ± 0.32 
-3.45 ± 0.54 


2.7/2 (25.9%) 
2.0/4 (73.6 %) 
2.19/1 (13.9%) 



Table 5.3: The parameters obtained for power law spectra (see Eq. \5.3\) for different phase intervals. 



lower energy observations. Most of the excess events are concentrated on SIZE < 100 whereas 
some excess events are still visible at SIZE > 100. The energy spectra have been calculated 
with a few different methods and all of them have given consistent results. All of PI, P2 and 
PI + P2 can be described by a power law from 25 GeV to 100 GeV and P2 has twice as high a 
flux as PI at 30 GeV. The power law indices of PI and P2 are compatible and they are approx- 
imately —3.5. Significant time variation of the pulsation between Cycles III and IV is not seen. 
These results will be discussed in more detail in Chapter [7] by comparing with Fermi-LAT data 
from 100 MeV to ~ 30 GeV. 
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Figure 5.20: The energy spectra calculated with the various methods. Red and green lines show the power 
law spectra obtained by forward unfolding based on SIZE distribution and the reconstructed energy 
distribution, respectively. Black, blue and pink points are the spectral points obtained by the unfolding 
with Tikhonov, Schmelling and Bertero regularization methods, respectively. Black lines indicate the 
power law fitting to the Tikhonov results while the shaded areas indicate the error of the fitting. 



Chapter 6 

Analysis of Fermi-L,AT public Data 



MAGIC could observe gamma-rays from the Crab pulsar only above 25 GeV, which is apparently 
beyond the spectral cut-off point. In order to make progress in understanding the emission mech- 
anism, the MAGIC results need to be discussed in connection with the lower energies. However, 
even after 9 years (April 1991 to May 2000) of operation of EGRET, which was the only GeV 
gamma-ray detector that could detect the Crab pulsar before 2007, only ~ 20 photons above 5 
GeV had been detected from the Crab pulsar (see Fig. ??). This gap in energy coverage between 
MAGIC and EGRET was soon to be filled by data from a new instrument. 

The Large Area Telescope (LAT) aboard the Fermi Gamma-ray Space Telescope was suc- 
cessfully launched on June 11, 2008. It can observe gamma-rays above 100 MeV and clearly 
saw the Crab pulsar up to ~30 GeV after 8 months of operation (see [4]). It is certain that the 
data of Fermi-LPiX help to interpret the results of MAGIC observations discussed in the previous 
chapter. Therefore, I analyzed one year of its data, which was made public in August 2009. 

6.1 Detector Design of Fermi-LAT 

The Fermi Gamma-ray Space Telescope is equipped with the Gamma-ray Burst Monitor (GBM) 
and the Large Area Telescope (LAT) (see the left panel of Fig. 16.11) . The GBM consists of 
12 thallium-activated sodium iodide (Nal(Tl)) scintillation counters and two bismuth germanate 
(BGO) scintillation counters (see H132II ). Each counter has an area of 126 cm 2 and energy ranges 
are 8 keV to 1 MeV and 200 keV to 40 MeV, for the Nal(Tl) counter and the BGO counter, re- 
spectively. The primary aim of the GBM is to detect gamma-ray bursts and its data have not been 
made public. The LAT comprises trackers, calorimeters, and an anti-coincidence detector (see 
OH). The LAT estimates the incoming direction and the energy of a gamma-ray by converting 
it into an electron-positron pair, which subsequently cause electromagnetic cascades inside the 
detector. 

Tracker: There are 16 tracker modules (see the right panel of Fig. 16.1) and each of them has 
18 (x, y) tracking planes. A tracking plane consists of two layers (x and y) of single-sided 
silicon strip detectors (35 cm long). The 16 planes at the top of the tracker are interleaved 
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Tracker 




coincidence 
lector 



■Calorimeter 



Figure 6.1: Left: A photograph of the Fermi Gamma-ray Space Telescope, adopted from HI 32]] . The high 
energy gamma-ray (> 100 MeV) detector LAT is seen at the top of the photograph. Six of the Nal(Tl) 
scintillation counters (for 8 IceV to 1 MeV) and one of the BGO scintillation counters (for 200 keV to 40 
MeV) which are for the gamma-ray burst monitor can also be seen. Right: The LAT detector. It consists 
of the tracker, the calorimeter and the anti-coincidence detector. Figure adopted from H31\l . See text for 
details. 



with high-Z converter material (tungsten) in which gamma-rays can convert to an electron- 
positron pair (see fl3TJ). The total vertical depth of the tracker including support material 
amounts to 1.5 radiation lengths. 

Calorimeter: The primary purpose of the calorimeter is twofold: 1) to measure the energy de- 
position due to the electromagnetic particle shower that results from the electron-positron 
pair produced by the incident photon and 2) to image the shower development profile, 
thereby providing an important background discriminator and an estimator of the shower 
energy leakage fluctuations. There are 16 calorimeter modules under the trackers (see the 
right panel of Fig. 16.11) . Each of them has 96 cesium iodide crystals doped with thallium 
(CsI(Tl)) with a size of 2.7 cm x 2.0 cm x 32.6 cm. The crystals are optically isolated 
from each other and are arranged horizontally in eight layers of 12 crystals each. The total 
vertical depth of the calorimeter is 8.6 radiation lengths, i.e. 10.1 radiation lengths for the 
total instrument (see PD ). 

Anti-coincidence detector: The purpose of the anti-coincidence detector is to provide charged- 
particle background rejection, u It surrounds the tracker modules (see the right panel of 
Fig. 16.11) . High energy gamma-rays may cause a so-called "backsplash" effect in the 



'The efficiency of charged particle detection is > 0.9997 (0.99999 when combined with the other subsystems). 
This is one of the key advantages compared to ground-based instruments, where it is impossible to install a primary 
hadron veto. 



6.2 Detector Performance of Fermi-LAT 
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Parameter Value or Range 

Energy range 20 MeV - 300 GeV 

Effective Area at normal incidence 9,500 cm 2 
Energy resolution (equivalent Gaussian la) 

100MeV-lGeV(on-axis) 9%-15% 

1 GeV-lOGeV(on-axis) 8%-9% 

10 GeV-300GeV (on-axis) 8.5%-18% 

> 10 GeV (> 60° incidence) < 6% 
Single Photon Angular resolution 

on-axis, 68% containment radius 6> 68 %: 

> 10 GeV <0.15° 
1 GeV 0.6° 
100 MeV 3.5° 

on-axis, 95% containment radius: 3 x0 68 % 

off-axis, containment radius at 55° 1.7 x on-axis value 

Field of View (FoV) 2.4 sr 

Timing accuracy < 10/is 

Event read-out time (dead time) 26.5 /is 

Table 6.1: Summary of LAT Instrument Parameters and Estimated Performance. 



massive calorimeter, i.e. isotropically distributed secondary particles (mostly 100-1000 
keV photons) from the electromagnetic cascade can hit the anti-coincidence detector (the 
photons cause Compton scattering), creating false veto signals. The effect was present in 
EGRET and lowered the detection efficiency above 10 GeV by at least a factor of two. To 
minimize the false veto, the anti-coincidence detector is made up of 89 segmented plastic 
scintillators: A 5 x 5 array on the top and 4x4 arrays on each of the four sides. Each 
tile is read out by two photomultipliers coupled to wavelength shifting fibers. Tiles near 
the incident candidate photon trajectory may be considered for background rejection (see 

mm). 



6.2 Detector Performance of Fermi-L AT 

The basic performance of the LAT is summarized in Table I6TI Parameters are taken from fl3~Tl . 
The energy resolution is dependent on the energy and the incident angle (see the table) but it 
is roughly 10%, which is better than MAGIC (35% at 30 GeV, see Fig. ??). The large FoV 
(2.4 sr) would also be of great advantage to compensate for the small effective area. The timing 
accuracy of < 10 //s is also good enough to study the Crab pulsar. On the other hand, the angular 
resolution is rather poor (0.6 degree at 1 GeV), leading to large contamination of the galactic 
diffuse emission or nearby source emissions to the target source. 

The systematic uncertainty in the energy scale was conservatively estimated to be < 5% for 
100 MeV to 1 GeV and < 7% above 1 GeV, from the comparison between electron beam tests 
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and their simulation (see [O and [|2)). The systematic uncertainty in effective area was evaluated 
by comparing Vela observation results and the simulation for them (see [|3)). It is 10% below 100 
MeV, decreasing to 5% at 560 MeV and increasing to 20% at 10 GeV and above. 

6.3 Data Sample 

After one year of operation, all the Fermi-LAT data and its analysis tools were made public in 
August 2009. I analyzed one year of data from 4th August 2008 to 3rd August 2009. Events 
with an energy between 100 MeV to 300 GeV and with an arrival direction of 20 degrees around 
the Crab pulsar were downloaded from the public Fermi website H210II . In order to have solid 
results, the events with a zenith angle smaller than 105 degrees and with the highest quality 
"Diffuse class", which means a high probability of being a photon, were selected. This event 
selection was performed by the Fermi official tool gtselect (see H21010 . Events with imperfect 
spacecraft information and events taken when the satellite was in the South Atlantic Anomaly 
were also rejected. These selections were carried out by gtmktime (see H21010 . 



Photon Map | Light Curve 




Right Ascension [deg] Phase 

Figure 6.2: Left: An event distribution of all the Fermi-LAT data used in the analysis (100 MeV to 300 
GeV) in equatorial coordinates. A star, a cross and a triangle indicate the positions of the Crab pulsar, IC 
443 and the Geminga pulsar. Emissions from these sources are clearly seen. In addition, galactic diffuse 
emission is also visible. Right: A light curve of the Crab pulsar based on the Fermi-LAT data from 100 
MeV to 300 GeV. The data between phase 0.52 and 0.87 (shaded area) are used for the nebula analysis. 
The nebula emission is a major background component for the pulsar analysis. 

The sky map of all the data used is shown in the left panel of Fig. 16.21 The Crab pulsar is at 
the center of the map (star) and two more bright sources, namely, IC 443 (cross) and Geminga 
(triangle), can also be seen. In addition, the galactic diffuse gamma-ray emission is visible. Since 
the angular resolution of the LAT detector is ~ 3.5 degree at 100 MeV (see PTlO . contamination 
from the nearby sources and the galactic diffused emission must be taken into account when the 



6.4 The Light Curve 
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Crab pulsar is analyzed. Moreover, emission from the Crab nebula must also be subtracted by 
using the pulse phase information. 

6.4 The Light Curve 

To make a light curve (phase diagram) of the Crab pulsar, first of all, the pulse phase must be 
assigned to each individual event. This is done by the official Fermi analysis tool, gtpphase (see 
[210J). It requires the pulsar ephemeris information and I used the "Jodrell Bank Crab Pulsar 
Monthly Ephemeris" for that, as I did for the MAGIC analysis (see Sect ??). Some of phase- 
resolved sky maps are shown in Fig. ??. Then, events from the direction around the Crab pulsar 
were extracted from the data set. Because the angular resolution has a strong energy-dependency, 
the extraction radius R [degree] should also be energy-dependent. I used the following radius R, 
which was also used in the Fermi official publication (see flU): 

R = Max(6.68 - 1.76 * log 10 (£), 1.3)) (6.1) 

where E is the estimated energy which is already assigned for each event in the public data. R 
decreases linearly to log 10 (£') until 1.14 GeV and stays constant at 1.3 degrees above this energy. 
It should be noted that the emissions from both the Crab pulsar and the Crab nebula are included 
within R. The angular resolution of Fermz-LAT does not enable a spatial resolution for the pulsar 
and the nebula. 

The light curve obtained by all the LAT data used is shown in the right panel of Fig. 16.21 
Energy-dependent light curves from 100 MeV to above 10 GeV are shown in Fig. 16.31 Below 10 
GeV, the pulsations are seen with good precision. A flat continuum in the light curves is mainly 
from the continuous Crab nebula emission. Above 10 GeV, although the two pulses are clearly 
visible, the statistical uncertainties are larger than for lower energies. PI (phases -0.06 to 0.04) 
and P2 (phases 0.32 to 0.43) have 12.1 ± 6.4 and 20.9 ± 7.3 excess events above 10 GeV with 
a significance of 2.0 and 3.2, respectively, with the background level (mainly from the nebula 
emission) estimated using the phases between 0.52 and 0.87. 

There are several features visible in these light curves: The flux ratio between PI and P2 is 
changing with the energy.The widths of the two pulses are decreasing as the energy goes higher. 
A hint of a possible third peak is visible at a phase around 0.75, but only above 10 GeV. These 
detailed features of the light curves will be discussed in Chapter [7] together with the MAGIC 
results and lower energy observations. 
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Figure 6.3: The light curves for different energies. From the top: 100 - 300 MeV, 300 - 1000 MeV, 1-3 
GeV, 3-10 GeV and above 10 GeV. The energy dependence ofPl/P2 ratio and the pulse width are visible. 
The possible third peak is also seen at a phase around 0.75, only above 10 GeV. Aflat continuum in the 
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6.5 The Crab Nebula Analysis 

As is the case with MAGIC, the rather poor angular resolution of Fermz'-LAT does not enable 
a spatial resolution for the Crab pulsar and the Crab nebula. Therefore, the nebula emission is 
the dominant background of the Crab pulsar emission for Fermi-LAtT data. It must be properly 
analyzed by using the pulse phase information and, then, must be subtracted from the Crab pulsar 
emission. Fortunately, the analysis of the nebula has an important by-product: one can make sure 
of the analysis method by checking if the obtained energy spectrum is smoothly connected to the 
IACT measurements (see Sect. ??). 

In order to analyze the nebula component, photons with the pulse phases from 0.52 to 0.1 
(see the right panel of Fig. 16.21 ), where no pulsed emission is seen in lower energies, are assumed 
to be from the nebula (and other background photons such as galactic diffuse emission). The 
selection of the right phase events is carried out by the Fermi-hAT analysis tool gtselect. The 
effective observation time and the collection area are calculated by gtltcube and gtexpmap (see 
[210]). The spectrum is determined by means of the likelihood method, using the official tool 
gtlike (see H210II ). It is done in the following way: The spectral shapes with several parameters of 
the sources in the FoV, the galactic and extragalactic diffuse emission models and the detector re- 
sponse function are assumed a priori. Then, the best parameters that maximize the likelihood of 
the observed data sets are determined. P6_V3_Diffuse, which is officially provided by the Fermi- 
LAT collaboration, was used for the detector response function, A simple power law spectrum 
was assumed for the IC 443 while a power law spectrum with an exponential cut-off was as- 
sumed for Geminga. For the extragalactic and galactic diffuse emission, isotropicJem_v02.txt 
and glLiem_v02.fit, which are included in the Fermz'-LAT analysis tool package as a standard 
model, were used. For the Crab nebula, the spectrum based on the sum of the two power laws is 
assumed aiming for the synchrotron and the inverse Compton emission components, which have 
been suggested by the previous EGRET (see H115II ) and IACT measurements (see e.g. 11201 . ifTTll 
and Fig. ??). 

The spectrum of the Crab nebula calculated based on the Fermz'-LAT data is shown in Fig. 16.41 
as a red line. One can see that synchrotron spectrum is steeply falling from 200 MeV to 500 MeV 
and the inverse Compton component becomes dominant above 800 MeV. They can be described 
as 

- di ^ 1 = (8.6 ±1.4) x lO^OE/SOOMeVr 4 - 1 ^ - 47 (6.2) 

dEdAdT sync 

(7.3 ± 0.7 x 1(T 12 (F/1 GeVr L68±a05 (6.3) 



dEdAdT ic 

(6.4) 



In order to make sure that the assumption of the spectral shape of the Crab pulsar is valid, 
the same data sets were divided into many subsets according to the energy. Then, the likelihood 
analysis was applied to each subset, assuming a simple power law in each small energy range. 

2 It is almost the same as the OP (off -pulse) phases defined in Sect. ?? but not exactly. Since the official Fermi 
publication used these phases for the nebula analysis (see PI). I followed their example. 
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The red points in the figure indicate the results for the divided subsets. Instead of showing many 
short truncated lines, the value at the bin center and its error are shown. All the points are very 
well aligned along the line, showing the validity of the assumption. 




Figure 6.4: The energy spectrum of the Crab nebula. The red line and circles indicate the spectrum 
calculated by myself using one year of Fermi-LAT data. The green dashed line indicates the published 
spectrum from the Fermi collaboration ( see /@/J. Blue circles, black open squares and black filled squares 
indicate the spectrum calculated by myself with the MAGIC observation data with the SUM trigger, the 
published spectrum from the MAGIC collaboration before the SUM trigger was installed ( see [20 / ) and the 
published spectrum from the HESS collaboration (see fliil/ ). respectively. The synchrotron and the inverse 
Compton components below 100 GeV are indicated by black dotted and black dashed lines, respectively. 



The spectrum published by the Fermi-LAT collaboration with smaller data samples ( 8 months 
of data) (see [4]) are also shown in the same figure as a green dashed line, which is consistent 
with my analysis. The measurements by MAGIC (the published one and the one I calculated 
with the data samples used for the pulsar analysis) and by HESS are also shown in the same 
figure. The spectra are smoothly connected from 100 MeV to above 10 TeV.0 

It should be noted that the poor statistics of Fermi-LAT data in the overlapped energy region 
from 50 GeV to 300 GeV does not allow a relative flux scale calibration between Fermi-LAT and 
MAGIC with a precision better than ~ 60%, which is larger than the systematic uncertainties of 
both experiments (see Sect. ?? and Sect. 16.21) . 

3 The true spectrum should not be a simple power law from 10 to 100 GeV. Therefore, the ~ 50% difference 
between the fitting line of Fermi-LAT data and that of MAGIC data at 100 GeV in Fig. |6.4| is not very meaningful. 



6.6 The Crab Pulsar Spectrum 
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6.6 The Crab Pulsar Spectrum 



From the nebula analysis, all the background components (the nearby source emissions, the 
diffuse emissions and the nebula emission itself) have been determined. Keeping the spectral 
parameters for background components fixed, the pulsed component is analyzed with the same 
likelihood method. Analyses are made for four different phase intervals; TP (total pulse, phase 
0.00 to 1.00), PI (-0.06 to 0.04), P2 (0.32 to 0.43) and PI + P2 (sum of PI and P2). These 
definitions are the same as those described in Sect. ??. 

6.6.1 Power Law + Cut-off Assumption for the Crab Pulsar Spectrum 

Previous EGRET measurements show that the energy spectrum of the Crab pulsar can be well 
described by a power law between 100 MeV to a few GeV (see H115II ). whereas non-detection by 
IACTs above 100 GeV (see EOl ) imply the cut-off to be around 10 GeV as mentioned in Sect. 
??. Therefore, I assumed the spectral shape as 



There are 4 free parameters, namely, the flux at 1 GeV / x , the cut-off energy E c , the power law 
index Ti and the sharpness of the cut-off T 2 0- As discussed in Sect ??, in the case where the 
emission region is close to the pulsar surface, the cut-off should be sharp, i.e. T 2 should be sig- 
nificantly larger than 1, while, in the case where emission comes from the outer magnetosphere, 
T 2 should be 1 . Therefore, estimation of T 2 is important for the investigation of the emission 
mechanism. However, due to the lack of statistics above 10 GeV, the likelihood analysis with 
four free parameters sometimes gives unstable results. Therefore, I made five analyses with three 
free parameters, i.e. with ft, E c , and Ti being free parameters and with T 2 fixed to be 0.66, 1, 
1.33, 1.66 and 2. 

The results are shown in Fig. 16.51 Fig. I6.6l and Table 16 .21 Hereafter, when T 2 = 1, the spectral 
cut-off shape will be called the "exponential cut-off", while when T 2 = 0.66, it will be called the 
"sub-exponential cut-off". The rest (with T 2 > 1) will be called the "super-exponential cut-off". 
In order to evaluate which assumption is more appropriate, the likelihood ratio — \og(L/L ex ) is 
calculated and shown in the seventh column of Table 16.21 where L and L ex are the likelihood 
value for a given assumption and that for the exponential cut-off assumption, respectively. The 
eighth column shows the corresponding probability. The spectral parameters published by Fermi- 
LAT collaboration for TP under the exponential cut-off assumption are also shown in the last 
row of Table 16 .2] They are consistent with my analysis (compare with the first row). 

What one can see from Fig. 16.51 Fig. I6.6l and Table 16 .21 is the following: 

4 As can be seen in Eq. 16.51 although the flux above the cut-off energy is suppressed, a low flux can still be 
expected above the cut-off energy. A word "cut-off does not mean extinction of a gamma-ray flux. 

5 ferm/-LAT collaboration reported these spectral parameters only under the exponential cut-off assumption in 
their publication J4). In addition to the TP spectrum, they analyzed the spectrum of many narrow phase intervals 
(width ~ 0.01). However, they did not publish the spectrum of PI or P2, which I want to compare with MAGIC 
results in Chapter 7. 
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• No significant difference is seen in the power law index r\ among different phase intervals 
(although r\ is dependent on the sharpness of the cut-off T 2 due to a mathematical effect). 
For the exponential cut-off assumption, Ti ~ —2.0. 

• The super-exponential assumptions are significantly worse than the exponential cut-off 
one. On the other hand, the sub-exponential cut-off assumption is as good as the exponen- 
tial cut-off one. 

• 80% of the total flux is from PI + P2 at 1 GeV. PI has twice as high a flux as P2 at 1 GeV. 

• The cut-off energy is higher for P2 than for PI. The difference in the cut-off energy be- 
tween TP and PI + P2 implies a higher cut-off energy for the bridge emission than for PI 
+ P2. 

• The flux of PI and P2 become comparable at around 5 GeV because of the higher cut-off 
of P2. P2 dominates above 5 GeV. 

The fact that super-exponential assumptions lead to worse fitting than the exponential one 
suggests that the emission region of gamma-rays is not close to the neutron star surface. This 
will be discussed further in Chapter[7] For the exponential cut-off model, the cut-off energies are 
estimated to be 6.1 ± 0.5, 4.5 ± 0.3, 3.7 ± 0.3 and 5.9 ± 0.7 GeV for TP, PI + P2, PI and P2 
respectively. 
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Crab Pulsar Spectrum: Exponential Cutoff 




Energy [MeV] 



Figure 6.5: The energy spectra of the Crab pulsar for different phase intervals. The thick solid lines 
are determined by the likelihood analysis, assuming the spectral shape to be a power law with an ex- 
ponential cut-off. The points are obtained by applying the same likelihood analysis to the small energy 
range, assuming a power law energy spectrum within the range (see text for details). For the the first two 
points below 200 MeV, where the systematic error strongly dominates the statistical error, the quadratical 
convolution of the systematic and statistical errors are indicated as pink lines. 
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Figure 6.6: The same as Fig. \6.5\ but with different spectral shape assumptions, namely, sub-exponential 
cut-off (top left), super-exponential cut-off with T2 = 1.33 (top right), super-exponential cut-off with 
T2 = 1-66 (bottom left), and super-exponential cut-off with T2 = 2.00 (bottom right) 
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a Likelihood ratio defined as LR = — \og(L/L ex ), where L and L ex are the likelihood value for a given assumption and 

that for the exponential cut-off assumption. 
b The probability corresponding to the LR value. It is also expressed as the corresponding deviation in the Gaussian distri- 



bution. When LR is negative, the probability is not calculated. 

Table 6.2: The results of the likelihood analyses for different spectral shape assumptions. For the defini- 
tion of the spectral parameters, see Eq. 16.51 The corresponding spectra are graphically shown in Fig. 16.51 
and Fig. \6.6\ The last row shows the spectral parameters from the Fermi-LAT publication for TP phase 
under the exponential cut-off assumption, taken from For other cut-off assumptions and for other 
phases (PI, P2 or PI + P2), the spectral parameters were not reported in ^jj. 
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6.6.2 Power Law Extension Assumption for the Crab Pulsar Spectrum 
above the Cut-Off 

It should be noted that the sharpness of the cut-off T 2 is already determined below ~ 7 GeV 
and the spectral shape above 7 GeV is not well determined due to the lack of statistics, as can 
be seen from the error bars in Fig. 16.5 l and Fig. 16.61 Considering that the MAGIC results show 
power law spectra with an index of ~ —3.5 above 25 GeV (see Sect. 15.101) . I examine a power 
law assumption for the spectra above the cut-off energy. The fact that the last point at 20 GeV 
is upwardly deviated from the exponential cut-off spectrum by ~ la for all phase intervals (see 
Fig. 16.51) may support this assumption. 

I selected data above 4 GeV and made a spectral analysis assuming a power law; 

;MM=/lo ( £ / 1 0GeV) r (6.6) 

The results are shown in Fig. I6.7l and Table 1531 In order to perform the likelihood ratio test 
with respect to the exponential cut-off assumption, the exponential cut-off assumption is also 
applied to the same data set. This time Ti is fixed to the best value obtained by the previous 
analyses (see Table 16.21) , so that the number of free parameters are the same (two) for the two 
models. The likelihood ratio (LR) and the corresponding probability are shown in the fifth and 
sixth columns of the table. None of the LRs are significantly large, which means that above 
4 GeV the power law assumption is as good as the exponential cut-off assumption. The obtained 
spectral indices are ~ —3.3 ±0.2 and consistent with the ones obtained from the MAGIC data, as 
described in Sect. 15 .101 (see Table 1531 ) . The spectra obtained with MAGIC data and Fermi-h AT 
data will be compared in more detail in the next chapter. 
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a Likelihood ratio defined as LR = — \og(L/L ex ), where L and L ex are the likelihood value for the power law 

extension assumption and that for the exponential cut-off assumption. 
b The probability corresponding to the LR value. It is also expressed as the corresponding deviation in the 



Gaussian distribution. When LR is negative, the probability is not calculated. 

Table 6.3: The results of the likelihood analyses assuming a power law spectrum ( see Eq. \6.6\ for different 
phase intervals. The corresponding spectra are graphically shown in Fig. \6.7\ 



6.7 Concluding Remarks 

By using one year of the public Fermi-LKT data, the Crab pulsar was analyzed from 100 MeV 
to ~ 30 GeV. Energy-dependent light curves show a few remarkable features, such as energy 
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Crab Pulsar Spectrum: Power Law 




Energy [MeV] 



Figure 6.7: The energy spectra of the Crab pulsar above 4 GeVfor different phase intervals. The thick 
solid lines are determined by the likelihood analysis assuming a power law spectrum. The points are 
obtained by applying the same likelihood analysis to the small energy range. Above 4 GeV, the spectra 
can be well described by a power law. 

dependence of the pulse width and the P2/P1 ratio. The possible third peak is also seen only 
above 10 GeV. These features will be discussed in detail, combined with the MAGIC results 
and lower energy observations, in Chapter [7] The energy spectrum of the Crab pulsar can be 
well described by a power law with an exponential cut-off, indicating that the emission region 
is not close to the neutron star surface. The cut-off energies are estimated to be 6.1 ± 0.5 GeV, 
4.5 ± 0.3 GeV, 3.7 ± 0.3 GeV, and 5.9 ± 0.7 GeV, for the TP, PI + P2, PI and P2, respectively. 
Due to the lack of statistics, the spectral shape above 7 GeV is not well determined. If data only 
above 4 GeV are analyzed, a power law function with an index of ~ —3.3 ± 0.2 fits the data 
too, suggesting the possibility that the energy spectrum is extending by a power law after 
the cut-off. Actually, the spectral indices are consistent with the MAGIC results. These spectra 
obtained with Fermi-LPR data will be compared with the ones obtained by MAGIC in Chapter 

m 
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Chapter 7 



Analysis of the Energy Spectrum and the 
Light Curve Combining MAGIC and 
Fermi-L,AT results 



MAGIC could measure the energy spectrum of the Crab pulsar from 25 GeV to 100 GeV, whereas 
Fermi-LAI could measure from 100 MeV to ~ 30 GeV. In this chapter, the energy spectrum and 
the light curve of the Crab pulsar from 100 MeV to 100 GeV is further investigated, by combining 
the results of the two experiments. 

7.1 Energy Spectra of PI + P2, PI and P2 

The energy spectra measured by MAGIC and Fermi-LAI are shown in Fig. 17.11 The top, the 
bottom left and the bottom right panels show those for PI + P2, PI and P2, respectively (see 
Sect. ?? for the definition of PI and P2). Although the Fermi-LAI results are well described by 
a power law with an exponential cut-off or with a sub-exponential cut-off, the MAGIC results 
are apparently deviated from them. In Sect. 17.1.11 and Sect. 17 .1.21 these deviations will be 
quantitatively examined, taking into account the systematic uncertainties of both experiments. 

On the other hand, a power law can well describe both Fermi-LAI results above 4 GeV 
(see Sect. 16X21) and MAGIC results between 25 GeV and 100 GeV (see Sect. I5T01 The 
obtained spectral slopes from the two experiments seem similar. The power law assumption 
will be examined in Sect. 17.1.31 also taking into account the systematic uncertainties of both 
experiments. 



7. Analysis of the Energy Spectrum and the Light Curve Combining MAGIC and 
54 Fermi-LAT results 




Figure 7.1: The comparison of the energy spectra of the Crab pulsar between Fenni-LAT measurements 
and MAGIC measurements for P 1 + P2 (top), PI (bottom left) and P2 (bottom right). The black and the 
red points indicate Fermi-LAT and MAGIC results, respectively. The sky blue and the green lines indicate 
the best fit spectra based on Fermi-LAT data, assuming a power law with an exponential cut-off and with a 
sub-exponential cut-off, respectively. The blue and the pink lines show the power law fits to the Fermi-LAT 
data above 4 GeV and MAGIC data, respectively. The statistical uncertainties of the fits are also shown 
by butterfly shape boxes with the corresponding color. Power law fits to the Fermi-LAT-MAGIC combined 
data above 4 GeV are shown by black dotted lines (see Sect. \7.1.3i . 
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7.1.1 How Much Do the MAGIC Measurements Deviate from an Expo- 
nential Cut-off Spectrum? 

In general, the theoretical calculations based on the outer magnetosphere emission hypothesis 
(the OG model or the SG model) predict an exponential cut-off at 1 - 10 GeV in the Crab pulsar 
energy spectrum, as discussed in Sect. ??, with which Fermi-LAT results are indeed consistent. 
However, as can be seen from Fig. 17. 1[ the MAGIC measurements seem to be deviated from 
the extrapolation from the Fermi-LAT measurements assuming that the spectra of PI, P2, and 
PI + P2 above 25 GeV follow a power law with an exponential cut-off. Here, I evaluate these 
deviations. 

Method: \ 2 Test on SIZE Distributions 

The power law spectra with an exponential cut-off for PI, P2 and PI + P2 obtained by Fermi- 
LAT data (see Table 16.21) are assumed to be also valid above 25 GeV. By means of MC with 
these energy spectra, expected SIZE distributions of excess events for PI, P2 and PI +P2 in the 
MAGIC data are computed. Then, x 2 tests are performed between these MC predictions and the 
actual observed distributions. The binning of the SIZE distribution is 30 — 50, 50 — 100 and 
100 — 400, which is the same as for the spectrum calculations (see Sect. 15.10.11) except that the 
highest two bins are combined so that all bins have a meaningful number of excess events (with 
respect to the statistical uncertainty) in the data. Among the free parameters in the Fermi-LAT 
likelihood analysis, the uncertainty of the cut-off energy has the largest effect on the expected 
SIZE distribution in the MAGIC data. Therefore, the analysis is repeated while changing the 
cut-off energy from 1 GeV to 25 GeV. This method is schematically shown in Fig. 17.21 

It should be noted that this method does not involve the energy reconstruction, which suffers 
the poor energy resolution and the trigger bias effect below 100 GeV (see Sect.??). SIZE (total 
charge in a shower image) is one of the most reliable image parameters and a good indicator of 
the gamma-ray energy. Therefore, a robust evaluation of the deviation can be expected for this 
method. 
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Figure 7.2: Schematical explanation of the method for evaluating the inconsistency between MAGIC and 
Fermi-LAT under the exponential cut-off spectrum assumption. A power law with an exponential cut-off is 
assumed for the spectral shape in MC. A power law part is taken from the Fermi-LAT measurements, while 
the cut-off energy is varied from 1 GeV to 25 GeV. A x 2 test is performed between the SIZE distribution 
of the observed data and that of the MCfor each cut-off energy. See text for details. 



Results of the x 2 Tests on SIZE Distributions 

The top left panel of Fig. 17.31 shows the x 2 value as a function of the cut-off energy. The 
number of degree of freedom is three. The corresponding upper probability is expressed with 
the corresponding Gaussian deviation and is indicated by black dotted lines. The green lines on 
the plot indicate the cut-off energies with statistical errors obtained from the Fermi-LAT data, 
i.e. (4.45 ± 0.31) GeV , (3.68 ± 0.29) GeV and (5.86 ± 0.74) GeV for PI + P2, PI and P2, 
respectively (see Table 16.21) . The top right panel of Fig. 17.31 shows the comparison of SIZE 
distributions between the observed data and the MC predictions for the Fermz'-LAT-determined 
cut-off energies. Based on these SIZE distributions, the inconsistencies between the exponential 
cut-off spectra determined by Fermi-LAT and the MAGIC measurements are estimated to be 
at (6.77 ± 0.13)(7, (3.01 ± 0.06)a and (6.04 ± 0.26)a levels for PI + P2, PI and P2, 
respectively. The cut-off energies that minimize the x 2 values are estimated to be (11.7 ± 0.7) 
GeV, (8.9 ± 0.9) GeV and (15.4 ± 1.2) GeV for PI + P2, PI and P2, respectively, which are also 
in clear contradiction with the Fermz'-LAT-determined cut-off energies. The bottom panel shows 
the comparison of SIZE distributions between the observed data and the MC predictions for the 
cut-off energies that minimize the x 2 values. These results are summarized in Table 177X1 
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X 2 vs Cut-off Energy (Exponential Cut-off) 
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Figure 7.3: The evaluation of the consistency between the Fermi-LAT and the MAGIC results under the 
exponential cut-off assumption. Top left: y 2 value based on the SIZE distribution of the MC prediction 
and that of the observed data, as a function of the cut-off energy. The number of degree of freedom is 
three. The corresponding upper probability is expressed with the corresponding Gaussian deviation and 
is indicated by black dotted lines. The black, the red and the blue points indicate PI + P2, PI and P2, 
respectively. Green bars on the plot indicate the cut-off energy and its statistical uncertainty obtained with 
the Fermi-LAT data. Top Right: The comparison of the SIZE distributions between the observed data 
and MC predictions. The cut-off energies obtained from the Fermi-LAT data are used for MC. Bottom: 
The same as the top right panel but the cut-off energies which minimize the x 2 ar e used for MC. 
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Effect of Systematic Uncertainties 

As mentioned in Sect. ?? and I6.2[ the systematic uncertainties in the energy scale of MAGIC 
and Fermi-LKT above 10 GeV are 20% and 7%, respectively, while the systematic uncertainties 
in the effective area (including the effective observation time) of MAGIC and Fermi-LKT above 
10 GeV are 10% and 20%, respectively. The apparent contradiction between Fermi-LAT and 
MAGIC measurements under the exponential cut-off assumption could be not due to the wrong 
assumption but due to these systematic uncertainties of the two experiments. 

In order to examine this possibility, I made the same analysis but scaling down the energy for 
MAGIC by 30%, i.e., the spectrum F MAGIC (E) used in the MC is defined as 

F MAGIC (E) = F LAT (0.7E) (7.1) 

where F LAT (E) is the power law with an exponential cut-off whose power law part is based 
on the Fermi-LAT measurements and whose cut-off energy is a parameter ranging from 1 GeV 
to 25 GeV. This 30% would be very conservative compared to the systematic uncertainties of 
both experiments. On the other hand, no correction for the effective area is applied because the 
uncertainty of the energy scale should dominate the effect on the results, due to the steep fall-off 
of the spectrum at MAGIC energies. 

The x 2 value as a function of the cut-off energy is shown in the left panel of Fig. 17.41 The 
comparison of the SIZE distributions between the observed data and the MC predictions for the 
Fermz'-LAT-determined cut-off energies is shown in the right panel of the figure. Even after 30% 
of energy scaling, the discrepancies in the SIZE distributions are (5.63±0.33)o\ (2.32±0.25)<r 
and (4.86 ± 0.60)cr level for PI + P2, PI and P2, respectively.The cut-off energies that minimize 
the x 2 values are (8.2 ± 0.5) GeV, (6.1 ± 0.6) GeV and (10.9 ± 0.9) GeV. These values are also 
significantly inconsistent with the Fermz'-LAT-determined values. These results are summarized 
in Table O 

From these discrepancies, it is evident that the exponential cut-off assumption is not valid 
at MAGIC energies (> 25 GeV). The extension of the pulsed gamma-ray emission observed by 
MAGIC requires a modification of the current pulsar models, which will be further discussed in 
Sect. IO 



7.1 Energy Spectra of PI + P2, PI and P2 



59 



X 2 vs Cut-off Energy (Exponential Cut-off) 
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Figure 7.4: Left: x 2 value based on the SIZE distribution of the MC prediction and that of the observed 
data, as a function of the cut-off energy. The number of degree of freedom is three. The corresponding 
upper probability is expressed with the corresponding Gaussian deviation and is indicated by black dotted 
lines. The black, the red and the blue points indicate PI + P2, PI and P2, respectively. Green bars on 
the plot indicate the cut-off energy and its statistical uncertainty obtained with the Fermi-LAT data. The 
energy scale of the MAGIC data are artificially lowered by 30% in order to examine the possibility that 
systematic uncertainties of both experiments is the reason for the discrepancy. Right: The comparison 
of the SIZE distributions between the observed data and the MC predictions for the cut-off energies 
obtained by Fermi-LAT data. 
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7.1.2 How Much Do the MAGIC Measurements Deviate from an Sub- 
Exponential Cut-off Spectrum? 

As discussed in Sect.?? and Sect. ??, nearly monoenergetic electrons (as many of theoretical 
models assume, see e.g. 119311 , H178H . H105H and H176H ) produce an exponential cut-off, i.e., T 2 in 
Eq. I6.5l is 1 . On the other hand, if T 2 is smaller than 1 , the inconsistency between the Fermi-LKT 
measurements and the MAGIC measurements would become smaller than for the exponential 
cut-off assumption. 

Since Fermi-LAT results can also be explained by the sub-exponential cut-off assumption 
(r 2 = 0.66, see Sect. 16.6.11 ), I examine the sub-exponential cut-off assumption. It should be 
noted that T 2 < 1 implies that the responsible electrons are not monoenergetic. In this case, no 
theoretical predictions for T 2 exist. I take the sub-exponential cut-off assumption (T 2 = 0.66) as 
one of the most extreme assumptions Q 

Method: x 2 Test on SIZE Distributions 

The method is the same as the one for the exponential cut-off assumption (see Sect. 17.1.1b except 
that T 2 = 0.66 (see Eq. |6.5b is used for the steepness of the cut-off. 

Results of the x 2 Tests on SIZE Distributions 

The top left panel of Fig. 17.31 shows the x 2 value as a function of the cut-off energy. The 
number of degree of freedom is three. The corresponding upper probability is expressed with the 
corresponding Gaussian deviation and is indicated by black dotted lines. The green lines on the 
plot indicate the cut-off energies for the sub-exponential cut-off with statistical errors obtained 
from the Fermi-LAT data, i.e. (2.20 ± 0.20) GeV , (1.79 ± 0.20) GeV and (3.00 ± 0.51) GeV 
for PI + P2, PI and P2, respectively (see Table l6\2l . The top right panel of Fig. 17.51 shows the 
comparison of SIZE distributions between the observed data and the MC predictions for the 
Fermz'-LAT-determined cut-off energies. Although discrepancies are smaller than those for the 
exponential cut-off case because of the more gradual cut-off, the inconsistencies between the 
sub-exponential cut-off spectra determined by Fermi-LAT and the MAGIC measurements are at 
(6.04 ± 0.26)cr, (2.55 ± 0.23)a and (5.35 ± 0.58)<r levels for PI + P2, PI and P2, respectively. 
The cut-off energies that minimize the x 2 values are (4.8 ± 0.3) GeV, (3.4 ± 0.4) GeV and 
(6.6 ± 0.5) GeV for PI + P2, PI and P2, respectively, which are also in clear contradiction 
with the Fermz-LAT-determined cut-off energies. The bottom panel shows the comparison of 
SIZE distributions between the observed data and the MC predictions for the cut-off energies 
that minimize the x 2 values. 

'in the case of the Vela pulsar, the energy spectrum of the total pulse (phase to 1) measured by Fermi-LAT 
shows T2 in Eq. 16 . 5 1 to be 0.69 ± 0.02 (see 15)). On the other hand, the energy spectrum from each narrow phase 
interval (width ~ 0.01) is still consistent with an exponential cut-off spectrum. In this analysis of the Crab pulsar, 
the widths of phase intervals are 0.1, 0.1 1 and 0.21 for PI, P2 and PI + P2. Therefore, T2 = 0.66 would be rather 
an extreme assumption. 
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Effect of the systematic uncertainties 

The effect of the systematic uncertainties of the two experiments can be taken into account in the 
same way as for the exponential cut-off case, i.e. the spectrum Fmagic{E) used in the MC is 
scaled down by 30% (see Eq. 17.11) . 

The x 2 value as a function of the cut-off energy is shown in the left panel of Fig. 17.61 The 
comparison of the SIZE distributions between the observed data and the MC predictions for the 
Fermz'-LAT-determined cut-off energies are shown in the right panel of Fig. 17.61 After 30% of 
energy scaling, the inconsistencies in the SIZE distribution are at (4.37±0.65)<r, (1.25±0.62)a 
and (3.71 ± 1.08)<r levels for PI + P2, PI and P2, respectively. The sub-exponential cut-off 
spectrum might explain the energy spectrum of PI if the systematic uncertainties of the two 
experiments are conservatively taken into account. However, PI + P2 and P2 are still largely 
inconsistent. The cut-off energies that minimize the \ 2 values are (3.39±0.23) GeV, (2.40±0.26) 
GeV and (4.77 ± 0.38) GeV. These values are still inconsistent with the Fermz-LAT-determined 
values. These results are summarized in Table I7TI 

The power law with a sub-exponential cut-off is not valid for PI + P2 and P2 at the MAGIC 
energies (> 25 GeV). 
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Figure 7.5: The same as Fig. \7.3\ but for the sub-exponential cut-off assumption. 
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Figure 7.6: The same as Fig. \7.4\ but for the sub-exponential cut-off assumption 
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Model 


Phase 


E Fermi rQ e V] a 


E* 1 [GeV] b 


Scale c 


X 2 from SIZE d 


Rejection Power e 


Exponential 
Cut-off 


PI +P2 


4.45 ±0.31 


11.68±0.74 




53.4-57.0 


6.6-6.9(7 


PI 


3.68 ± 0.29 


8.86± 0.91 




15.3-16.4 


3.0-3.1 a 


P2 


5.86 ±0.74 


15.39± 1.19 




42.1 -48.8 


5.8 -6.3 a 


Sub- 
exponential 
Cut-off 


PI +P2 


2.20 ± 0.20 


4.77±0.31 




42.1 -50.2 


5.8-6.4a 


PI 


1.79 ±0.20 


3.43±0.36 




11.3-14.2 


2.3-2.8a 


P2 


3.00 ±0.51 


6.59± 0.52 




30.8-44.0 


4.8-5.9a 




Exponential 
Cut-off 


PI +P2 


4.45 ±0.31 


8.17±0.53 


0.7 


36.5-44.4 


5.3-6.0(7 


PI 


3.68 ± 0.29 


6.12± 0.63 


0.7 


9.9-12.8 


2.1-2.6(7 


P2 


5.86 ±0.74 


10.86± 0.85 


0.7 


25.8-38.3 


4.3-5.5(7 


Sub- 
exponential 
Cut-off 


PI +P2 


2.20 ± 0.20 


3.39± 0.23 


0.7 


21.1-33.5 


3.7-5.0(7 


PI 


1.79 ±0.20 


2.40± 0.26 


0.7 


4.0-8.9 


0.6-1.9(7 


P2 


3.00 ±0.51 


4.77±0.38 


0.7 


13.2-31.0 


2.6-4.8a 



a The cut-off energy determined by the Fermi-LPd data 

b The cut-off energy estimated by the SIZE distribution in MAGIC data, assuming that the power law part determined by 

the Fermi-LAI data is valid. 
c The energy scaling factor. This factor is applied to the MAGIC energy. 

d x 2 value calculated from the SIZE distributions of the observed data and the MC prediction. The number of degree of 
freedom is three. 

e The probability that the two results are consistent, expressed in the corresponding Gaussian deviation. This is calculated 
from the \ 2 value. 

Table 7.1: Evaluation of the inconsistency between the MAGIC measurements and the exponential/sub- 
exponential cut-off spectra determined by Fermi-LAT. 
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7.1.3 Power Law Fit above 4 GeV 

The super-exponential cut-off assumption is ruled out by the Fermz'-LAT results alone (see Sect. 
16.6.11) . Moreover, it is now evident that neither the exponential cut-off nor the sub-exponential 
cut-off can explain the observational results of Fermi-LAI and MAGIC consistently. 

On the other hand, a power law can well describe both the Fermi-LAT results above 4 GeV 
(see Sect. I6A21) and MAGIC results between 25 GeV and 100 GeV (see Sect. 15101 as can be 
seen in Fig. 17.11 In this section, an assumption that the spectra above 4 GeV follow a power law 
is examined. 



Method: Combined Fit to Fermi-LAT and MAGIC Data 

A power law function is fitted to the data points of Fermi-LAT and MAGIC together. In order to 
take into account a possible energy scale difference between the two experiments, the correction 
factor k is introduced to the power law function as: 

a = f 1 (E < 25GeV) 

\ k (E> 25GeV) K } 

Fermi-LAT points are statistically independent of each other but MAGIC data points are 
correlated, due to the unfolding procedure. This correlation will be taken into account when 
fitting is applied. Firstly, k is fixed to 1 and then, secondly, k is treated as a free parameter. 



Results of Combined Fit to Fermi-LAT and MAGIC Data 

The results are shown in Fig. 17 .7 l and Table 17 .21 Even without relative energy correction, i.e. with 
k fixed to 1, the fitting probability is high enough (see the fifth column of the Table). By making 
K a free parameter, %2 values are only slightly reduced and the obtained ks do not significantly 
deviate from 1 . As a result, f w and T do not change significantly. T is about 3 for all the phase 
intervals and a significant difference is not seen. 

The physics interpretation of this result will be discussed in Sect J8.21 



Phase 


/io [10~ 7 cm-^TeV" 1 ] 


r 


Scaling Factor k 


X 2 /dof(Prob.) 


PI +P2 


2.96 ± 0.20 


-2.96 ± 0.08 


1 


8.09/5(0.15) 


PI 


1.26 ±0.16 


-3.18 ±0.16 


1 


1.98/5 (0.85) 


P2 


1.63 ±0.13 


-2.81 ±0.08 


1 


9.29/5(0.10) 


PI +P2 


2.61 ±0.30 


-3.12 ±0.14 


0.83 ± 0.09 


5.54 / 4 (0.24) 


PI 


1.16 ±0.21 


-3.29 ± 0.22 


0.87 ±0.15 


1.40/4(0.84) 


P2 


1.40 ±0.21 


-3.00 ±0.18 


0.82 ±0.12 


7.56/4(0.11) 



Table 7.2: The power law fit combining the Fermi-LAT data above 4 GeV and the MAGIC data 
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Figure 7. 7: The power law fitting, combining the Fermi-LAT data above 4 GeV and the MAGIC data. The 
correlation among the MAGIC points is taken into account. The left panels show the case when the energy 
scales of both experiments are used while the right panels show the case when a shift in the energy scale 
of MAGIC is allowed by letting k in Eq. \7.2\ be a free parameter. Top, middle and bottom panels show 
PI + P2, PI and P2, respectively. The best fit parameters and fitting probabilities are shown in Table \7.2] 
The fitted lines without the energy scaling (left panels) are also shown in Fig. \7.I\ as black dotted lines. 
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Figure 7.8: Overlaid light curves from different energies. The vertical values are normalized such that 
the heights of the PI peak are equal. The energy dependence of the P2/P1 ratio and Bridge/Pi ratio are 
clearly seen. 



Name 


PI 


Bridge 


P2 


OP 


Phase Interval 


-0.96 to 0.04 


0.14 to 0.25 


0.32 to 0.43 


0.52 to 0.88 



Table 7.3: The Definition of Names of Pulse Phases 



The light curves for different energies such as Fig. 15 . 1 5 1 and Fig. 16.31 suggest that the flux 
ratio between PI and P2 changes with energy. The fraction of the Bridge emission seems to be 
energy-dependent, too. Some of the light curves from different energies are overlaid in Fig. 17.81 
in order to show the energy dependence of the P2/P1 ratio and the Bridge/Pi ratio. The P2/P1 
ratio and/or the Bridge/Pi ratio have been studied by many authors such as H115II . H184II . 111331 . 
Hi 271 and H128B for a wide energy range from optical to gamma- rays. On the one hand, the 
energy dependence of these ratios can be thought of as a consequence of the different energy 
spectra for different phase intervals. On the other hand, taking the flux ratio between two phase 
intervals cancels out the systematic uncertainty of the absolute flux scale for different detectors. 
Therefore, these ratios enable to study precisely the relative spectral behavior of different phase 
intervals for a very wide energy range using the measurements from many different detectors. 

Here, I calculate the P2/P1 and Bridge/Pi ratios for Fermi-LKT and MAGIC data and com- 
pare them with the lower energy results. 
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Method: Calculation of Energy-dependent P2/P1 and Bridge/Pi Ratios from Light Curves 

Fermi-LAT data sets are divided into 5 bins in energy, 100 - 300 MeV, 300 - 1000 MeV, 1 - 3 
GeV, 3 - 10 GeV and above 10 GeV. MAGIC data are divided into two bins in SIZE, 25 - 70 
ph.e. and 70 - 500 ph.e.. The number of excess events for PI, P2 and Bridge are calculated from 
light curves of each energy bin, by estimating background by OP phases (For the definition of 
PI, P2, Bridge and OP, see Sect. ?? and Table O)- Then, the ratio of P2 to PI and that of 
Bridge to PI are computed. 

Results of Energy-dependent P2/P1 and Bridge/Pi Ratios 

The numbers of excess events of PI, P2 and Bridge as a function of energy are shown in the 
left panel of Fig. 17.91 and Table [7~4l In the energy ranges of above 10 GeV for Fermi-LAT and 
70 - 500 ph.e. for MAGIC, no significant excess was found in Bridge phases. Therefore, 95% 
confidence level upper limits are shown for these bins. The horizontal central value of each point 
corresponds to the log-mean energy of each bin, taking into account the detector's effective area 
and the energy spectrum of the pulsar, while the horizontal error bar corresponds to 30% of the 
log-mean energy, which is a conservative estimation for the absolute energy uncertainty. The 
P2/P1 ratio stays almost constant from 100 MeV to 3GeV and rises at energies above 3 GeV. On 
the other hand, the Bridge/Pi ratio rises from 100 MeV to 3 GeV. Above 3 GeV, due to the lack 
of statistics, it is not possible to draw a conclusion in the behavior. 

In Fig. 17. 10L the results are compared with those at lower energies. The data points of lower 
wavelengths are adopted from H115I . A few common features between the behavior of the P2/P1 
ratio and that of the Bridge/Pi ratio can be seen: From 1 eV to 1 MeV, they increase with energy 
with a power law. From 1 MeV to 100 MeV, they drop rapidly and, then, rise again above 
100 MeV - 1 GeV. The ratios at ~ 30 GeV may be as high as those at 1 MeV, although the 
uncertainty is large due to the lack of statistics. 



Energy 


Energy a 


PI 


P2 


Bridge 


P2/P1 


Bright/Pi 


Range 


[GeV] 


[counts] 


[counts] 


[counts] 






Fermi-LAT 














0.1 -0.3 GeV 


0.156 


4135 ±71 


2067 ± 55 


103 ± 33 


0.50 ± 0.02 


0.025 ± 0.008 


0.3-1.0 GeV 


0.481 


5219 ±76 


2885 ± 59 


254 ± 28 


0.55 ± 0.01 


0.049 ± 0.005 


1.0-3.0 GeV 


1.52 


1706 ± 42 


974 ± 33 


224 ± 18 


0.57 ± 0.02 


0.131 ±0.011 


3.0-10.0 GeV 


4.38 


248 ± 17 


245 ± 17 


79 ± 11 


0.99 ±0.10 


0.318 ±0.050 


> 10 GeV 


12.7 


12 ±6 


21 ±7 


< 17 b 


1.70 ± 1.05 


< 1.34 


MAGIC 














25 - 70 ph.e 


24.0 


4711 ± 1129 


8233 ± 1198 


1338 ± 1195 


1.75 ±0.49 


0.284 ± 0.261 


70 - 500 ph.e. 


51.2 


1437 ± 899 


3096 ± 954 


< 2953 b 


2.16 ± 1.50 


< 2.06 



a The representative energy, taking into account the pulsar spectrum and the detector response. 
b 95% confidence level upper limit. 



Table 7.4: The numbers of excess events of PI, P2 and Bridge and their ratios 
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Figure 7.9: Left: The number of excess events as a function of energy. The red, blue, and green points 
denote PI, P2 and Bridge, respectively. The last two points above 20 GeV are based on MAGIC data, 
while the rest are based on Fermi-LAT data. Right: The P2/P1 ratio (blue) and the Bridge/Pi ratio 
(green), as a function of energy. The last two points above 20 GeV are based on MAGIC data, while the 
rest are based on Fermi-LAT data 
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Figure 7.10: The P2/P1 ratio (left) and the Bridge/Pi ratio (right) for a wide energy range from 1 eV to 
100 GeV. Black points are adopted from HI 15V . Green points and a green arrow are based on Fermi-LAT 
data, while red points and a red arrow are based on MAGIC data 
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7.3 Rising and Falling Edges 

The energy spectra have been calculated for the specific phase intervals. They are especially 
important to deduce the emission regions of the pulsar on a large scale. On the other hand, as 
can be seen from Fig. I5.15l and Fig. 16.31 the pulse shape is not the same for all energies. A study 
of the change of pulse shape with energy would also be very helpful to understand the emission 
mechanism in more detail, i.e. on a small scale, because the pulse phase should be connected to 
the geometry of the emission region, as discussed in Sect. ??. 

Due to the lack of statistics (with respect to the large background), it is difficult to analyze 
the pulse shape with MAGIC data alone. However, the fact that the pulse peak phases are very 
similar for all energies implies that the pulsations for all energies may originate from a common 
physical process. Therefore, one can also expect common features in pulse shapes. Once a 
common feature is found, one can statistically examine it in the MAGIC data. 

In fact, I found that by plotting the light curves in log scale, some interesting features become 
visible. In Fig. 17.1 11 light curves for different energies from optical to very high energy gamma- 
rays are shown. The optical light curve is obtained from the MAGIC central pixel data. L. 
Kuiper provided the X-ray light curves, which have also been used in 111 1511 . They are based on 
ROSAT HRI (100 eV to 2.4 keV, see EH), BeppoSAX MECS (2.4 keV to 10 keV, see 112551 ). 
BeppoSAX PDS (20 keV to 100 keV, see 12051 ) and CGRO COMPTEL (750 keV to 30 MeV, see 
H2081I ). The gamma-ray light curves from 100 MeV to 10 GeV are produced by myself based on 
the Fermi-LAT data. The gamma-ray light curve at 20 to 200 GeV □ are taken from my analysis 
of the MAGIC data with SIZE between 25 to 500 ph.e.. One can see the following features: 

• Both rising and falling edges show an exponential behavior. 

• Slopes are not symmetric between rising and falling edges. 

• Slopes change with energy. 

Here, I discuss the pulse edges of PI and P2. In the previous sections and chapters, PI and 
P2 referred to specific phase intervals, PI being from -0.06 to 0.04 and P2 being from 0.32 to 
0.43, while in this section they simply denote the first peak and the second peak. 

Method: Fitting Exponential Functions to Pulse Edges 

An exponential function 

F(p) = A x exp(±(p - po)/r) (7.3) 

is fitted to rising and falling edges of PI and P2 for different energies, A, p and r being free 
parameters. Fitting ranges were chosen such that the bridge emission and the pulse peak do not 
worsen the goodness of the fit. Fitting ranges are summarized in Table 17.51 Then, the energy 
dependence of the rise time r r i se and the fall time Tf a ii for PI and P2 are examined. 



2 This energy range is a rough estimation based on the SIZE range from 25 to 500 ph.e.. No significant excess 
is seen at 200 GeV. 
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Figure 7.11: Light curves in a logarithmic scale for different energies. From the top left downward to 
the bottom left, optical measurements by the MAGIC central pixel (see Sect. 15.61 ), 100 eV to 2.4 keV by 
ROSAT-HR1, 2.4 keV to 10 keV by BeppoSAX MECS, 20 keV to 100 keV by BeppoSAX PDS, 0. 75 MeV to 
WMeVby COMPTEL and 10 MeV to 30 MeV by COMPTEL. From the top right downward to the bottom 
right, 100 MeV to 300 MeV, 300 MeV to 1 GeV, 1 GeV to 3 GeV, 3 GeV to 10 GeV measured by Fermi-LAT 
and 20 GeV to 200 GeV ( 25 < SIZE < 500) measured by MAGIC. 
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Edge 


Rise of PI 


Fall of PI 


Rise of P2 


Fall of P2 


Phase 


-0.07 to -0.005 


0.005 to 0.04 


0.29 to 0.395 


0.405 to 0.50 



Table 7.5: Fitting range for the rise/fall time estimation 



Results of the Exponential Function Fitting to Pulse Edges 

Fitted lines are overlaid on the light curves shown in Fig. 17. Ill The obtained parameters and 
the fitting probabilities are shown in Table 17.61 Fitting probabilities are reasonably good for 
all energies. The energy dependence of the rise and the fall time of both peaks are graphically 
shown in the upper panels of Fig. 17.121 The horizontal values in Fig. 17 .121 are determined in the 
following way: Below 100 MeV, horizontal error bars show the full energy range of the data sets 
and the central values are the logarithmic center of the range. Above 100 MeV, the horizontal 
central values and error bars are determined in the same way as the P2/P1 ratio study (see Sect. 

m. 

General behavior of the rise time and the fall time from optical to gamma-rays can be sum- 
marized as follows: 

• Below 10 MeV, the rise time of both pulses and the fall time of PI are increasing, while 
the fall time of P2 is decreasing. 

• Above 100 MeV, the rise time of PI and the fall time of P2 are decreasing, while the rise 
time of P2 and the fall time of PI do not show clear energy dependence. 

Hereafter, I would like to focus on the behavior above 100 MeV, where I have personally 
analyzed all the data. In the lower panels of Fig. I7.12[ the energy dependence of the rise and the 
fall times above 100 MeV are shown. They are well described by a linear function of logio(-E') 
or a constant value, where E is the energy. By fitting a linear function to the rise time of PI and 
the fall time of P2 and by fitting a constant value to the fall time of PI and the rise time of P2, 
the following is obtained: 



2 - (9.4 ± 1.3) x 10~ 3 log 10 (E[GeV]) (7.4) 
2 (7.5) 

(7.6) 

2 - (9.6 ± 3.1) x 10~ 3 log 10 (E[GeV]) (7.7) 

X 2/dof are 2.2/3, 12.2/4, 4.6/4 and 2.4/3 for r^ e , T™ e , rf^and rg„ respectively. 
MAGIC points and Fermi-LAT points can be fitted consistently. It is also interesting that the 
rise time of PI and the fall time of P2 show very similar dependence on energy. The energy 
dependence of the pulse shape above 100 MeV is schematically illustrated in Fig. 17.131 As 
energy increases, the outer edges, i.e. the rising edge of PI and the falling edge of P2 get sharper 
while the inner edges, i.e. the falling edge of PI and the rising edge of P2 do not change. 
The physical interpretation of these results will be discussed in Sect. 18.61 



r r Z(E) = (2.02 ±0.08) x 10- 

T ™ e = (6.46 ± 0.24) x 1(T 

Tfa = (1.73 ±0.08) x 10- 

t™(E) = (2.42 ± 0.16) x 10- 
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a This energy range is a rough estimation based on the SIZE range in MAGIC data from 25 to 500. No significant 
excess is seen at 200 GeV. 



Table 7.6: Results of the rise and the fall time estimation for P 1 and P2 
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Figure 7.12: Energy dependence of the rise and the fall times. Top left: The rise time of PI (red) and P2 
(blue), as a function of energy. Top right: The fall time of PI (red) and P2 (blue), as a function of energy. 
Bottom left: The same as the top left panel but the energy range from 100 MeV to 100 GeV is zoomed. 
Bottom right: The same as the top right panel but the energy range from 100 MeV to 100 GeV is zoomed. 




Figure 7.13: An illustration of the energy dependence of the pulse edges above 100 MeV. As energy 
increases, the outer edges get sharper, while the inner edges do not change. 
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7.4 Peak Phases 

As described in Sect. ??, the pulse profiles of the Crab pulsar are very similar at all energies, 
while a closer look at the light curves reveals a slight energy dependence of the peak phase. For 
example, the first peak in the X-ray data (see Fig. ??) and that in optical data (see Fig. 15.1 II) 
precede the radio peak by ~ 0.01 in phase. Also, above 100 MeV, there seems to be a slight 
shift of the peak phase, which can be seen if the light curves based on the Fermz-LAT data are 
zoomed, as shown in Fig. 17.141 



Fermi 100 -300 MeV 

— Fermi 1 - 3 GeV 

— Fermi 10 -300 GeV 




Phase Phase 



Figure 7.14: A closer look at the peak phase of PI (left) and P2 (right). The black, red and green lines 
indicate the light curves for 100 MeV to 300 MeV, 1 to 3 GeV and above 10 GeV, based on Fernii-LAT 
data. 

Here, I examine the energy dependence of the peak phase above 100 MeV by using Fermi- 
LAT and MAGIC data. Apparently, the energy dependence of the peak phases is not very strong. 
The precision of ~ 0.003 in phase would be required to discuss the energy dependence. In 
order to fulfill this requirement with the obtained data, a sophisticated method called the "Kernel 
density method" is used. 

7.4.1 Kernel Density Method 

The Basic Concept 

If the statistical significance of the signal is large enough, the peak phase can be precisely deter- 
mined without assuming a specific pulse shape. However, the significance of the obtained signal 
is not high enough to determine the peak phase with a precision of ~ 0.003, especially for ener- 
gies above 10 GeV. By assuming the pulse shape a priori, the fitting of the assumed shape to the 
obtained data might improve the precision. In such a case, in order to achieve the best possible 
precision, the data should not be converted into a binned phase diagram but should be analyzed 




Fermi 1 00 - 300 MeV P2ZOOITI 

Fermi 1 - 3 GeV 
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on an event-by-event basis. An event-by-event maximum likelihood method could have been 
used for that but, in the case of MAGIC data, which is dominated by ~ 15 million background 
events, it requires too much computational power. 

The Kernel density method also uses event-by-event information but it does not require too 
much computational power. Moreover, the assumption of the pulse shape is not needed either, 
although a so-called "kernel estimator" must be chosen beforehand. The kernel density method is 
a well established statistical method for estimating the probability density function of a measured 
parameter, based on the observed data sample (see e.g. [1801 , H179M and ||36l ). The true pulse 
profile can be interpreted as a probability density function for the pulse phase of the signals and, 
hence, the kernel density method can be applied. 

The probability density f(p) as a function of phase p can be estimated as follows: 

m - vi K *<cTr) (7 - 8) 

where N, p i7 Kh(x) and h are the total number of events, the phase of zth event, a kernel estimator 
and the band width of the kernel estimator K h . The method is schematically illustrated in Fig. 
17.151 In a light curve (a phase histogram), the phase of a given event is smeared by a kernel 
density estimator. f(p) is the sum of the smeared curves of all the events. An example of the 
application of the method using the MAGIC data is shown in Fig. 17.161 








Phase 



Figure 7.15: Schematical explanation of the Kernel density method. In a phase histogram, each event is 
smeared by a kernel density estimator, as shown in the left panel. f(p) is the sum of the smeared curves 
of all the events, as indicated by a red line in the right panel. 



The Choice of the Kernel Estimator K h and its Width h 

As a kernel estimator Kh, a normal Gaussian is commonly used in many different applications, 

Kh {—) = 7S e (7 - 9) 

On the other hand, as shown in Fig. 17.171 a Lorentzian represents the pulse shape better. The 
Lorentzian kernel estimator is written as 

(nr) - ^tt^e <7 ' 10) 

(7.11) 
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| MAGIC 25 < Size < 500 | | MAGIC 25 < Size < 500: Gaussian with medium width" 
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Figure 7.16: An example of the kernel density method. Left: The original light curve from MAGIC data 
with SIZE between 25 and 500. Right: The obtained probability density function. The Gaussian kernel 
with h = 0.024 is used. 



One of these estimators should be used. 

Not only the shape of the estimator but also the width h must be properly chosen (see H187I ). 
Too narrow a width will make the density function f(p) too wiggly and produce many spurious 
features. On the other hand, too big a width will lead to too smooth a function that smears out 
all the structure. Since the pulse shape is not symmetric, as described in the previous section, the 
smearing may cause a shift in the peak phase. I chose h based on the light curve of Fermi-LAT 
data above 3 GeV, which is more or less the (logarithmic) center of the concerned energy range. 
First, a Gaussian and a Lorentzian are fitted to PI and P2 independently. Results are shown in 
Fig. 177171 The best fit hs are 0.016 ± 0.002 (Gaussian fit to PI), 0.031 ± 0.004 (Gaussian fit to 
P2), 0.011 ± 0.001 (Lorentzian fit to PI) and 0.024 ± 0.003 (Lorentzian fit to P2). The optimal 
h should be close to these values. 
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Figure 7.17: Fitting a Lorentzian and a Gaussian for a light curve of Fermi-LAT data above 3 GeV 



The final choice of the shape of K h (Gaussian or Lorentzian) and the width h was made such 
that the peak phase shift caused by the smearing effect is minimal. The phase shift is estimated 
as follows: Assuming a pulse shape to be 

ftrueip) - I exp(-p/T /o „) if p>0 ^ 
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f(p) is calculated by convoluting the pulse shape ft rue (p) with the kernel estimator K h (p) 



0.25 



0.25 



ftrue(p')K h (p - p')dp' 



(7.13) 



Then, the difference in the peak phase Ap peak between ft rU e(p) an d f(p) is estimated. Hereafter, 
App ea k is referred to as the "analytical phase shift". An example of the analytical phase shift is 
shown in Fig. 17.181 
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Figure 7.18: Left: An example of the convolution of a true pulse with a kernel estimator. T r i se = 0.03 and 
T faii = 0.02 are used for the true pulse shape ftrue(p)> while a Gaussian kernel K^ip) with h = 0.016 
is used for the kernel estimator. A black and a red line indicate the true pulse shape ftrue(p) an d the 
convoluted pulse shape f(p), respectively, while a green line indicates the kernel estimator Kh(p) which 
is scaled such that the peak height = 1. Right: The same as the left panel but phases from -0.02 to 0.02 
are zoomed. The analytical phase shift Ap pea/ t is indicated by a blue arrow. 



By substituting Eqs. 17.41 — 17.71 for r rise and Tf a ii in Eq. 17.121 I estimated the energy de- 
pendence of the analytical phase shift Ap pea k for different K^s, as shown in Fig. 17. 191 The 
Lorentzian kernel estimators with h = 0.006 for PI and h = 0.012 for P2, which are the half of 
the fitted values (see Fig. 17.171) have the smallest effect. Even smaller h might reduce the effect 
further. However, if h is too small compared to the fitted values, /(p) would, in turn, produce 
spurious structures, as mentioned before. Therefore, in the following analysis, I use Lorentzian 
for Kh with h = 0.006 for PI and h = 0.012 for P2. The residual Ap pea k will be subtracted from 
the obtained results. 



Statistical Uncertainty Estimation by the Bootstrapping Method 

One can estimate a probability density function f(p) from observed data by the kernel density 
method and, then, determine the peak position from it. However, it does not give the statistical 
uncertainty. 

The statistical uncertainty of the result can be estimated by the bootstrapping method which is 
well established and is used in many statistical treatments (see e.g. Il70l and Il63l0 . The procedure 
is explained as follows: Let N be the total number of observed events. One randomly chooses N 
events out of the observed N events. The same events can be chosen multiple times. Then, the 
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Figure 7.19: The energy and the estimator dependence of the analytical phase shift Ap pea kfor PI (left) 
and P2 (right) 

same kernel density method is applied to the chosen data sample. By repeating this procedure 
M times, one obtains M different f(p)s. The RMS of the peak phase distribution from them is 
used as the statistical uncertainty of the peak phase. In this analysis, M = 900 is used. 



7. Analysis of the Energy Spectrum and the Light Curve Combining MAGIC and 
80 Fermi-LAT results 



7.4.2 Example of the Kernel Density for Different Energies 



In Fig. 17 .201 and Fig. 17.2 1 L the original light curves and the resulting probability functions (kernel 
densities) are shown. The statistical uncertainty of the function estimated by the bootstrapping 
method are indicated by colors. 
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Figure 7.20: The original light curves (left) and the resulting probability density function (kernel density) 
by using Lorentzian kernel with h = 0.012 (right). This kernel density is used for the P2 peak deter- 
mination, while h = 0.006 is used for the PI peak. Black dotted lines indicate the background level. 
900 curves obtained by the bootstrapping method are overlaid in the left panels and colors indicate their 
density. From the top: 100 - 300 MeV, 0.3 - 1 GeV and 1-3 GeV are shown, all of which are based on 
Fermi-LAT data. 
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Figure 7.21: The same as Fis \7.20\ but for higher energies. From the top: 3-10 GeV from Fermi-LAT 
data, above 10 GeV from Fermi-LAT data, 20 - 40 GeV (25-70 in SIZE) from MAGIC data and 40 - 200 
GeV (70 - 500 in SIZE) from MAGIC data. These energy ranges for MAGIC data are rough estimation 
based on the SIZE range and no significant excess is seen at 200 GeV. Because of the relatively low 
significance of the signals, an (acceptably) wiggly structure is visible. By using larger h in the kernel 
estimator, the structure will become less pronounced, while it leads to a larger analytical bias. Black 
lines in the left panels indicate the resulting probability density functions, while white lines are their (la) 
uncertainty. 
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7.4.3 Results 



Results are shown in Fig. 17.221 Open and filled squares indicate the results before and after 
correcting the analytical phase shift Ap p , respectively. Horizontal values and errors are deter- 
mined in the same way as the P2/P1 ratio study (see Sect. 17.21) . Blue lines indicate the energy 
dependence of Ap p estimated in Sect. 17.4. H Csee Fig. I7.19I ). 

After the correction (filled squares), the energy dependence is clearly seen for PI but it is 
not the case for P2, due to the large uncertainties. This difference comes mainly from the pulse 
width. Since P2 has twice as large a width as PI, the precision of the peak phase becomes worse. 
The results after the correction (filled squares) are fitted by a linear function and 



Peakl(E) 
Peak2(E) 



(-3.8 ± 0.6) x 10~ 3 + (2.1 ± 0.9) x 10^ 3 log 10 (£[GeV]) (7.14) 
0.39 + (3.9 ± 1.6) x 10" 3 - (0.05 ± 2.8) x 10~ 3 log 10 (E[GeV]) (7.15) 



are obtained, x^/dof are 6.26/5 and 3.18/5, for Peakl and Peak2, respectively. 
The physical interpretation of the results will be discussed in Sect. 18.71 
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Figure 7.22: The energy dependence of the peak phase of PI (left) and P2 ( right). Open and filled squares 
denote the peak phase before and after correcting the analytical phase shift Ap pea f,. Green lines indicate 
the Appeal (see Fig. \7.19i . The peak phases after the correction as a function of energy are fitted by a 
linear function, as shown by the black lines. 
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7.5 Possible Existence of the Third Peak 

As mentioned in Sect. 16.41 a possible third peak is seen above 10 GeV in Fermi-LAT data. 
Judging from the light curve (see Fig. 16.31 ), the flux might be as high as PI and P2, although the 
large statistical uncertainty does not permit a solid argument. Here, I examine the existence of 
the possible third peak by using both Fermi-LAI data and MAGIC data. 



7.5.1 Definition of the ON and OFF Phases for the Third Peak P3 

The third peak has not been detected in other energies, except for specific frequencies in radio, 
where fourth and fifth peaks are also seen (see Sect. ??). Therefore, it is not possible to define 
the third peak phase interval (P3) a priori. I define P3 to be phases from 0.7 to 0.8. It is based 
on the observed result itself, which leads to the overestimation of the flux and the statistical 
significance, while no fine tuning of the bin-edges is carried out in order to reduce the effect. In 
all the previous analyses, background level had been estimated by OP phases 0.52 to 0.88 (0.52 
to 0.87 for the Fermi-LAI nebula analysis). For P3, phases from 0.5 to 0.65 and from 0.85 to 0.9 
are used for the background estimation. 



7.5.2 MAGIC and Fermi-LAT above 10 GeV 
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Figure 7.23: The light curves of the Fermi-LAT data above 10 GeV (top) and that of the MAGIC data 
with SIZE from 25 to 500. P3 phases and background estimation phases are indicated by red and black 
shadows, respectively. 
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The light curves of Fermi-LAT data above 10 GeV and MAGIC data with SIZE from 25 
to 500 are shown in Fig. 17.231 Fermi-LAT data show 21.5 ±6.8 excess events corresponding to 
3.5cr, while MAGIC data show 2670 ± 1563 excess events corresponding to 1.7a. More statistics 
are required in order to verify or refute the presence of the signal. 

Although the existence of the signal is not clear, I estimated the energy spectrum of P3 based 
on Fermi-LAT data. Based on MAGIC data, the differential flux upper limit was also calculated 
with a 95% confidence level . They are shown in Fig. 17.241 The estimation is done in the 
following way: Instead of using the likelihood method for Fermi-LAT data and the unfolding 
method for MAGIC data, the number of excess events (or the upper limit on the number of 
excess events) in each energy bin was simply divided by the effective area and the observation 
time. Therefore, the spill-over effect from the adjacent bins are not taken into account. The 
discrepancy between MAGIC upper limits and Fermi-LAT measurements may imply an upward 
bias in the Fermi-LAT data analysis, probably because the P3 phases are defined based on the 
observed data themselves. The possibility that there is some signal, whose excess was enhanced 
by the background fluctuation, cannot be excluded. The time variability of the P3 excess is also 
another explanation. 
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Figure 7.24: The energy spectrum of P3. Black points and red arrows indicate the measurements by 
Fermi-LAT and the upper limits based on MAGIC data, respectively. The flux is estimated by simply 
dividing the number of excess events by the observation time and the effective area. A green and a blue 
line indicate the exponential cut-off energy spectra of PI and P2, respectively. 



7.6 Concluding Remarks 

Neither an exponential cut-off nor a sub-exponential cut-off can explain the energy spectra mea- 
sured by MAGIC and Fermi-LAT consistently, even if the possible systematic uncertainties of 
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both experiments are taken into account. Assuming that the energy spectrum does not roll off 
exponentially but extends with a power law after the cut-off energy, they can be consistently 
explained with an index of ~ —3.0 ±0.1 (after the cut-off). The physics interpretation of these 
spectral features will be discussed in Sect 18.21 

The P2/P1 ratio stays constant between 100 MeV and 3 GeV and rises above 3 GeV. On the 
other hand, the Bridge/Pi ratio increases by a power law between 100 MeV and 3 GeV, while 
the behavior above 3 GeV cannot be analyzed, due to the lack of statistics. The edges of the two 
pulses have an exponential shape. The rise time of PI and the fall time of P2 have a clear energy 
dependence while the fall time of P2 and the rise time of P2 do not. The physics interpretation 
of this will be discussed in Sect. 18 .61 The peak phase of PI has a slight energy dependence 
while this is not clear for P2, due to the difficulty in determining the P2 peak phase. The physics 
discussion on this will also be presented in Sect. 18.71 The third peak seen in Fermz'-LAT data 
above 10 GeV with 3.5 a is not clearly visible in MAGIC data. In order to verify or refute the 
existence of the P3 signal, more statistics is required. 



7. Analysis of the Energy Spectrum and the Light Curve Combining MAGIC and 
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Chapter 8 
Discussion 



Several interesting features of the pulsed gamma-ray radiation from the Crab pulsar above 100 
MeV have been newly found from the MAGIC and Ferm/-LAT data analyses. Quite a few results 
are not in agreement with model predictions and extrapolations. Especially the following findings 
are remarkable: 

• None of the super-exponential cut-off, the exponential cut-off and the sub-exponential cut- 
off can describe the measured energy spectrum. 

• The measured energy spectrum extends at least up to 100 GeV. 

• A power law with an index of —3.0 ± 0.1 can well explain the measured energy spectrum 
between 4 GeV and 100 GeV. 

• The edges of the two pulses have a clear exponential behavior. The rise time of PI and the 
fall time of P2 have a clear energy dependence while the fall time of PI and the rise time 
of P2 do not. 

• The peak phase of PI has a small but significant energy dependence. 

In this chapter, I discuss the new constraints on the pulsar emission models and possible modi- 
fications of the standard model, based on these findings. In addition, the radiation efficiency in 
gamma-rays above 100 MeV is discussed. 



8.1 Constraints on the Emission Region 

As discussed in ??, there are mainly two approaches to infer the emission region, namely, by the 
steepness of the cut-off and by the highest energy of the observed photons. 
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8. Discussion 



8.1.1 By the Steepness of the Cut-off 

As described in ??, if the emission region is close to the star surface, a superexponential cut-off is 
expected. T 2 in Eq. 16.51 is typically 2 (see e.g. ||141| and Il65l0 . However, the analysis of Fermi- 
LAT data revealed that superexponential assumption (T 2 = 2) is less likely than exponential 
cut-off (r 2 = 1) by 4.8a, 7.7a, 5.0a and 4.3a, for TP, PI + P2, PI, and P2 respectively. This 
indicates that emission region is far from the star surface and that gamma-rays do not cause 
magnetic pair creation. In essense, the modest steepness of the spectral cutoff allow one to reject 
the PC model. 

8.1.2 By the Highest Energy of the Observed Photons 

More quantitative estimation of the emission region can be made by the highest energy of the 
observed photons. As discussed in Sect. ??, the highest energy of photons E max which can 
escape from a given height r can be estimated as 

/ r \ 7/2 D 

E max {r) ~ 40VP — -f MeV (8.1) 

where P, R , B , and B cr are the period of the pulsar in second, the radius of the neutron star, 
the magnetic field strength at the stellar surface and the critical magnetic field (B cr = 4.4 x 10 13 
G). 

MAGIC detected gamma-rays up to 100 GeV. For P = 33.6 ms, B = 3.8 x 10 12 G = 
0.086£> cr and E max = 100 GeV, one obtains r/R = 7.8. The height of the emission region must 
be more than 7.8 times the pulsar radius, which is too large for the PC model. 

8.2 Estimates of the Electron Spectrum and 

Constraints on the Acceleration Electric Field 

From the argument in the previous section, it is clear that the emission region is free from mag- 
netic pair creation process, i.e. the emission region should be in the outer magneto sphere. Since 
magnetic pair creation is a basic physics process, this conclusion is rather robust. 

Even if the emission region is assumed to be in the outer magnetosphere, the energy spectrum 
observed by MAGIC contradicts the most favored theoretical models, which predict an exponen- 
tial cut-off. Below, I will discuss possible modifications of the standard outer magnetosphere 
model. After briefly reviewing the basic equations for the electron energy spectrum within the 
pulsar magnetosphere and the curvature radiation spectrum from these electrons in Sect. 18.2.11 
I will deduce the electron spectrum based on the observed gamma-ray spectrum of PI + P2 in 
Sect. 18.2.21 Then, constraints on the strength of the acceleration electric field will be discussed 
in Sect. 18.2.31 Possibilities to explain the deviation of the observational results from the standard 
model by a imperfect dipole magnetic field will be briefly described in Sect. 18.2.41 
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8.2.1 Basic Equations for the Electron Energy and the Curvature Radia- 
tion Photon Energy 

As discussed in Sect. ??, the energy of the electron is set by the equilibrium between the gain 
in energy due to the acceleration electric field and the energy loss via the curvature radiation. 
The Lorentz factor of the electron T can be written as a function of the acceleration electric field 
strength E\\ and the magnetic field curvature R curv : 

/ * w \ x / 4 



1 = 2 ' 8 X 10 /m7 J / \ V Rcurv /1000 km (8.2) 

y(10 7 eV/cm) J v 

See Sect. ??, for the derivation of this equation, 

As described in Sect. ??, the curvature radiation spectrum from a single electron with an 
energy Tm e c 2 is written as 

7 mono 

e c = h 3 hu curv = h 3 — *g— (8.4) 

Z Z Z1X rL curv 

From Eq. 18.21 and Eq. 18.41 the cut-off energy E c of the curvature radiation is written as a 
function of En and R curv : 

eEv ^ 3/4 



E c = 6.5 777^— rf-, — r jRcurv/WWkm [GeV] (8.5) 
y (10' eV/cm) ) v 
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8. Discussion 



8.2.2 Estimates of the Electron Spectrum 

Based on the Measured Gamma-ray Spectrum 

Here, I deduce the electron spectrum (the T spectrum) based on the measured gamma-ray spec- 
trum of PI + P2, assuming that gamma-rays above 4 GeV are generated by the curvature radia- 
tion. For simplicity, the curvature of the magnetic field is assumed to be fixed at R curv = 1000 
kmQ. 

0) Standard Model: Nearly Monoenergetic Electron Spectrum 

In general, it is considered that En and R curv of the last closed field line do not change extremely 
largely over the emission region (see e.g. H175I ). In addition, the dependency of T on En and 
R curv is rather weak (see Eq. I8.21 ). These are the reasons why a nearly monoenergetic spectrum 
is derived for accelerarted electrons in most of pulsar models based on the outer magnetosphere 
emission hypothesis. 

The cut-off in the high energy gamma-ray spectrum is determined by the curvature radiation 
spectrum from these nearly monoenergetic electrons, in turn leading to the exponential cut-off 
(see Sect. 18.2.11) . Examples of the theoretical explanations of the Crab pulsar spectrum observed 
before 2007 are shown in Fig. 18.11 which are the same as Fig. ??. The highest end of the 
spectrum is explaned by the curvature radiation from nearly monoenergetic electrons. In order 
to explain the power law spectrum with an index of ~ 2 between ~ 10 MeV and ~ 1 GeV, the 
synchrotron radiation (left panel) or the inverse Compton scattering (right panel) is considered. 
Below 10 MeV, the emissions from secondary electron-positron pairs created by high energy 
gamma-rays are considered to be responsible for the observed spectrum. 

Applying the standard scenario to the observed spectrum of PI + P2 around the cut-off energy 
is shown in Fig. 18.21 I made a simple calculation assuming that the T spectrum of the electrons 
has a Gaussian shape with the mean of 2 x 10 7 and the RMS of 10 6 . This corresponds to 
En = 3 x 10 6 [V/cm], which is consistent with, e.g. H175H . The absolute flux scale was chosen 
such that the predicted gamma-ray spectrum matches with the measurements. 

As can be clearly seen in the figure, the measured spectrum requires modification of the 
standard models. 



'The different R curv for different emission region can be taken into account by replacing T with 

r(i? c „™/1000[km])- 1 / 3 (see Ea.WM 
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Figure 8.1: Theoretical explanations of the energy spectrum of the Crab pulsar assuming that the emission 
region is in the outer magnetosphere. Two models, namely Harding et al. (left, see H93\l ) and Tang, et al 
(right, see H178\l ) are shown. The highest end of the spectrum are explained by the curvature radiation 
from nearly monoenergetic electrons as indicated by a solid line (left) and a dashed-dotted line ( right). The 
emission between 100 MeV and 1 GeV, where a power law spectrum with an index of ^ —2.0 is observed, 
is explained by the synchrotron radiation (a dash-dotted line in the left panel) or inverse Compton emission 
(a dotted line in the right panel). Below 10 MeV, the emissions from secondary electron-positron pairs 
created by high energy gamma-rays are considered to be responsible for the observed spectrum. 



Electron Spectrum Gamma-ray Spectrum: P1 + P2 




Figure 8.2: Applying the standard scenario to the highest end of the observed spectrum for PI + P2. Left: 
Assumed V factor spectrum for electron. It has a Gaussian shape with the mean of2x 10 7 and the RMS of 
10 6 . Right: The energy spectrum of the curvature radiation from the assumed electron spectrum, overlaid 
with the observed results. MAGIC measurements are largely deviated from the expected spectrum. The 
discrepancy between the curvature radiation spectrum and observed data below 3 GeV (a green shaded 
region) could be explained by either the synchrotron radiation or the inverse Compton scattering (see. 
Fig. H3. 
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8. Discussion 



1) Steep Power Law Tail Assumption 

Given the fact that the power law with an index of ~ —3 describes well the observed gamma- 
ray spectrum above 4 GeV, the simplest assumption for the electron spectrum would be a power 
law with an index of ~ —8, because if the electron has a power law spectrum with an index 
of —p, the resulting curvature radiation spectrum should show the power law with an index of 
—q = —(p + l)/3, as discussed in Sect ??. I assumed a power law with an index of —8 between 
1.7 x 10 7 and 10 8 in the T spectrum. For T < 10 7 , a Gaussian with the mean of 1.7 x 10 7 and 
the RMS of 10 6 is assumed so that the modification with respect to the standard model is only 
on higher energy side. The measured spectrum of PI 2 can well be explained by the curvature 
radiation from these electrons as expected (the top panels of Fig. 18.3b . Obviously, the main 
difficulty is to find a convincing argument for the shape of the electron spectrum. 

2) Log-Gaussian Assumption 

It might also be possible that the observed power-law-like behavior above the cut-off energy is 
part of a curved spectrum, whose curvature is not visible due to the statistical uncertainty, the 
limited energy resolution and the limited energy coverage of the measurements. 

For example, a log-Gaussian spectrum of the electrons produce a gamma-ray spectrum ex- 
tending to the MAGIC energies. In the bottom panels of Fig. 18.31 the log-Gaussian spectrum of 
the electron with the mean of log 10 (T) mean = 7 and the RMS of log 10 (T) RMS = 0.15 is assumed. 



The spectrum above the cut-off can be reasonably explained, taking into account that MAGIC 
and Fermi-LPJ' may have a relative energy scale difference of up to ~ 30%. The idea behind the 
log-Gaussian spectrum of T is that the electrons are not as monoenergetic as the standard model 
predicts. This might originate from a small distortion of the pulsar magnetosphere structure from 
the standard model. 




(8.6) 
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Electron Spectrum Gamma-ray Spectrum: P1 + P2 




Electron Spectrum Gamma-ray Spectrum: P1 + 92 




Figure 8.3: The same as Fig. \8.2\ but for different F spectrum. Top: A power law tail after the peak 
is assumed for the F spectrum. At F < 10 7 , a Gaussian with the mean of 1.7 x 10 7 and the RMS of 
10 6 is assumed, while a power law spectrum with an index of —8 is assumed between 1.7 x 10 7 to 10 8 . 
The resulting curvature radiation reproduces the observed results very well. Bottom: A Log-Gaussian 
spectrum is assumed for the F spectrum. The mean oflogio(F) mean = 7 and the RMS oflogio(F)nMS = 
0.15 is used. The resulting curvature radiation reproduces the observed results reasonably well taking 
into account the possible energy scale difference between MAGIC and Fermi -LAT up to ~ 30%. 
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8. Discussion 



3) Two Population Assumption 

Several authors consider the possibility that the observed pulsation is the sum of the emissions 
from the two poles (see e.g. H175B and H178I - ). The two contributions to the light curve are shown 
in Fig. 18.41 Even though their calculations do not predict an energy spectrum extending to 100 
GeV, the possibility of the contributions from the two isolated places make the two population 
assumption intriguing. 

In fact, the model of two population of electrons can explain the measurements as well, as 
shown in the top panels of Fig. 18.51 Here I assumed two Gaussian spectra with the mean T being 
2 x 10 7 and 4 x 10 7 . R curv = 1000 km is used. The RMSs of the two Gaussians are 20 times 
smaller than its mean. The peak flux of the second population is 200 times lower than that of the 
first. 
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Figure 8.4: A theoretical explanation of the light curve assuming that emissions from the two poles con- 
tribute to the observed pulsation. The inclination angle of the dipole axis with respect to the rotation axis 
is assumed to be 50 degrees and the viewing angle is assumed to be 76 degrees. Contributions from each 
pole are overlaid with different intensity of the lines. Figure adopted from M78\l 



4) Power Law + Exponential Cut-off Assumption 

The power law with an exponential cut-off for the electron spectrum produces an interesting re- 
sults. By assuming a power law with an index of —3 and cut-off at T = 0.5 x 10 7 for electrons, 
the resulting curvature radiation can explain the observed gamma-ray spectrum very well from 
100 MeV to 100 GeV, if the relative energy scale difference of up to ~ 30% is taken into account 
(see the bottom panels of Fig. 18.51) . This assumption does not require additional emission mech- 
anisms such as the inverse Compton scattering and the synchrotron radiation in order to explain 
the measured spectrum between 100 MeV to 3 GeV. 

It is known that the electron energy spectrum may exhibit a power law with an exponential 
cut-off in the case of shock acceleration in a supernova remnant taking into account either the 
acceleration-time limits (see |]68l) or radiative-loss limits (see H193I0 . In order to apply a sim- 
ilar scheme to a pulsar, the acceleration mechanism in the pulsar magnetosphere needs to be 
reconsidered from scratch. 
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Electron Spectrum Gamma-ray Spectrum: P1 + P2 




Electron Spectrum Gamma-ray Spectrum: P1 + 92 




Figure 8.5: The same as Fig. \8.2\ but for different T spectrum. Top: Two population of the T spectrum is 
assumed. Two Gaussian spectra with the mean T being 2 x 10 7 and 4 x 10 7 are used. The RMSs of the two 
Gaussians are 20 times smaller than its mean. The peak flux of the second population is 200 times lower 
than that of the first. The resulting curvature radiation reproduces the observed results very well. Bottom: 
The power law with an exponential cut-off is assumed for V spectrum. The power law index of —3 and 
cut-off at r = 0.5 x 10 7 are used. The resulting curvature radiation can explain the observed spectrum 
very well from 100 MeV to 100 GeV, if the relative energy scale difference of up to ~ 30% is taken into 
account. This assumption does not require additional emission such as the inverse Compton scattering or 
the synchrotron radiation between 100 MeV to 3 GeV. 
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8. Discussion 



8.2.3 Constraints on the Acceleration Electric Field in an Ideal Dipole 
Magnetic Field 



'°9 10 (n 
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Figure 8.6: The electron Lorentz. factor T as a function of the field curvature R CU rv and the acceleration 
field strength En, expressed with a color scale. Black lines indicate the contor lines for T = 2.0 7 , 4.0 x 10 7 
and 7.0 x 10 . Blue lines indicate the corresponding cut-off energy of curvature radiation photons for 5 
GeV, 20 GeV and 100 GeV (see Eq. 18.41 ). A red dotted line indicate the co-rotation radius, which should 
be the upper limit in R CU rv 



Based on Eq. 18.21 T as a function of R cur v and E\\ is graphically shown in Fig. 18.61 R CU rv 
as a function of En for E c = 5 GeV, 20 GeV and 100 GeV are shown by blue lines in the same 
figure. 

It is known that R CU rv of the dipole magnetic field is well approximated as R CU rv(r) = \[R~iJ 
(see e.g. 113410 . where Rl is the co-rotation radius (see Sect. ??) and r is the distance from the 
center of the neutron star. Since the emission region must be within the light cylinder, R curv < 
Rl would be the reasonable upper limit on R curv . Assuming that the (true) gamma-ray spectrum 
is extending as a power law up to 100 GeV without a cutoff, as is the case for 1) Steep Power 
Law Tail Assumption in the previous subsection, there must be a place where E\\ > 2.9 x 
10 8 [V/cm]. This is ~ 100 times larger than the value used in the standard model (see e.g. H175I0 . 
Assuming that the gamma-ray spectrum observed at around ~ 100 GeV is basically the tail of 
the curvature radiation spectrum with a 20 GeV cut-off, as is the case for 3) Two Population 
Assumption in the previous subsection, the lower limit in E» will be Eu > 3.3 x 10 7 [V/cm]. 
This is ~ 10 times larger than the standard value. 
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8.2.4 Possible Explanations for the High Energy Tail of the Gamma-ray 
Spectrum in an Imperfect Dipole Magnetic Field 



Figure 8. 7: Schematical explanations of two ideas for the high energy tail of the observed gamma-ray 
spectrum. Left: The magnetic field curvature near the light cylinder can be much larger than that of the 
dipole fields, because some field lines in the gap should connect to the wind zone. Such a large field 
curvature may produce a high energy tail in the gamma-ray spectrum as can be understood from Eq. 18.51 
Right: If the magnetic field lines are slightly wiggling along the dipole structure, there may be two radially 
connected regions the inner one of which has a larger-than-dipole field curvature and the outer one of 
which has a smaller-than-dipole field curvature. Since the curvature radiation cooling is weaker in the 
larger curvature field, electrons can be accelerated to higher energies. Then, these electrons subsequently 
enter the smaller curvature field, generating the unusually high energy curvature photons. 

The high energy tail of the observed gamma-ray spectrum could be explained by the imper- 
fect dipole structure of the magnetic field. 

The thickness of the Outer Gap or the Slot Gap can be as much as 0.\R L at the light cylinder. 
The magnetic field lines in the middle of the gap should not close within the light cylinder and 
should be connected to the wind zone. Therefore, there is a possibility that the magnetic field 
curvature R cur v near the light cylinder in the middle of the gap is significantly larger than that 
of the dipole structure (see the left panel of Fig. 18 .7b . As can be seen from Eq. 18.21 the larger 
the curvature (R CU rv)i the higher the electron energy (r). Alghouth the energy of the curvature 
radiation photon is proportinal to R~^ rv , it is also proportional to T 3 (see Eq. 18.41) . resulting in 
the higher photon energy from the larger curvature (see Eq. 18.51) . 

A more efficient mechanism for producing the high energy tail could be as follows: Let us 
assume that the magnetic fields are slight wiggling along the dipole field in a small scale (1 — 10 
km). There may be two radially connected regions the inner one of which has a larger-than- 
dipole field curvature and the outer one of which has a smaller-than-dipole field curvature (see 
the right panel of Fig. 18.71) . Since the curvature radiation cooling is weaker, electrons can be 
accelerated to higher energies in the larger curvature field. These electrons subsequently enter 
the smaller curvature field, generating the anormally high energy curvature photons. 
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8. Discussion 



8.3 Inverse Compton Scattering as a Second Radiation Com- 
ponent 

As shown in the top panels Fig. 18.51 a second emission component can explain the deviation of 
the MAGIC measurements from the standard model. In the previous section, a second population 
of electrons is assumed. By introducing the inverse Compton scattering, a single monoenergetic 
electron population might also explain the MAGIC measurements. 

As discussed in Sect. ??, the energy of the radiated photon via the inverse Compton scattering 
is (see Eq. ??) 

£ 7 ~ er 2 (8.7) 

where e is the energy of the target photon and T is the Lorentz factor of the electron. In order to 
have a radiation peak at 30 GeV for an electron population with T = 2 x 10 7 , e should be peaked 
at 

e ~ EJY 2 ~ 10" 4 [eV] (8.8) 

For a thermal radiation, 10~ 4 eV corresponds to 1 K, which is by far lower than the stellar 
surface temperature (~ 10 6 K, see Sect ??). On the other hand, the CMB radiation corresponds to 
2.7 K in temperature, which is actually close to ~ 10~ 4 eV. However, it is known that the number 
density of CMB photons (400 cm~ 3 ) is much lower than the radio photons emitted within the 
pulsar magnetosphere (> 10 7 cm -3 ) (see H103I0 . 

Therefore, in order to reproduce the inverse Compton scattering spectrum peaked at around 
30 GeV, there must be a particular mechanism to produce a soft photon spectrum peaked at 
~ 1CT 4 eV, which is hard to imagine. Therefore, an effect of a simple inverse Compton scattering 
for the extended spectrum is highly unlikely. 



8.4 Radiation Efficiency 

The total energy loss of a pulsar, i.e. the spin down luminosity E can be estimated from the period 
and the time derivative of it, as described in Sect ??. For the Crab pulsar, it is E = 4.6 x 10 38 
erg/s. 

Here, I add an estimate of the radiation efficiency, which is a fraction of the spin down 
luminosity deposited in a given energy range. The estimation method is adopted from JU. When 
the observed flux is F obs , the luminosity L can be calculated as L = An fiiF obs D 2 where /n 
and D are the beaming angle factor and the distance from the pulsar to the Earth. f n depends 
on the inclination angle a of the magnetic dipole axis and the viewing angle £ (see Sect ??). 
For the Crab pulsar, assuming the SG model or the OG model, a ~ 70 degrees and £ ~ 60 
degrees well explain the light curves (see e.g. 1153) ). the nebula torus structure in X-ray (see 
COED), and the polarization of the optical pulsation (see H170II ). For a and £ near these values, 
fa ~ 1.0 is the good approximation according to [11921 . The distance is known to be D = 
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2.0 ± 0.2 kpc (see e.g. H185II ). From the measurements, F obs above 100 MeV is calculated to be 
(1.31 ± 0.03) x 10" 9 erg/s. 



Using E 
calculated as 



4.6 x 10 38 erg, D = 2.0 kpc and f n = 1.0, the efficiency r] above 100 MeV is 



V 



L 

E 



A7rf n F obs D' 
E 



1.36 ±0.04) x 10" 



(8.9) 



Only 0.13% of the total energy loss is due to the radiation above 100 MeV. As discussed in Sect. 
?? and Sect. ??, the rotation energy is carried away mostly by the pulsar wind. 

Fig. 18.81 shows rj(E) as a function of energy. Based on the combined analysis described in 
Sect. 17.1.31 X] for energies between 30 and 100 GeV is (8.4 ±1.1) x 10~ 6 . If the spectrum had 
shown the pure exponential cut-off, it would have been 5.li| 8 x 10~ 7 , which is a factor of ~ 10 
lower than the measurement. The discrepancy between the exponential cut-off and the MAGIC 
measurement corresponds to (8 x 10~ 6 ) /(1.4 x 10~ 3 ) = 6 x 10~ 3 = 0.6% of the radiation energy 
above 100 MeV. 
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Figure 8.8: The radiation efficiency r] for different energy intervals. Red points indicate the observed 
results. Below 10 GeV, TP spectrum with the exponential cut-off assumption is used while for the bin 
of 10 GeV to 30 GeV, a power law assumption above 4 GeV is used (see Sect. 16.61 ). For the bin of 30 
GeV to 100 GeV, combined analysis with an power law assumption for PI + P2 is used (see Sect. \7.1.3\) . 
The black point indicates the extrapolation of the exponential cut-off spectrum of TP based on Fermi-LAT 
measurement. 
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8. Discussion 



8.5 Remarks on the High Energy Tail of the Crab Pulsar En- 
ergy Spectrum 

Fermi-LAI measured the energy spectrum of the Crab pulsar with good precision up to ~ 10 
GeV, which is only a factor of ~ 2 higher than the exponential cut-off energy (E c in Eq. |6.5l) . The 
Fermz'-LAT-measured spectrum is consistent with the standard model. Therefore, the standard 
model is successful in explaining the energy spectrum of the Crab pulsar in the energy region 
where a vast majority of the gamma-ray radiation energy is deposited. 

On the other hand, MAGIC observed the energy spectrum of the Crab pulsar above 25 GeV, 
which is a factor of ~ 5 higher than the exponential cut-off energy. Only a small fraction of the 
radiation energy is deposited above 25 GeV. The deviation from the standard model detected by 
MAGIC requires only a higher order correction for the standard model. 

It should also be mentioned that, all the dedeuced electron spectra discussed in Sect. 18.2.21 
eventually need to be corrected by the "smearing" of the energy due to the limited and not- 
perfectly-known energy resolution. While the simulations clearly showed that the observed spec- 
trum up to 100 GeV is not just an artifact of the energy resolution (see Sect. 17.1.11) . I could not 
determine the precise partition of genuine high energy gamma-rays and lower energy gamma- 
rays mismeasured with higher energy assignment. It is obvious that clarification can only come 
from better measurements. 

8.6 Energy Dependence of the Rising and Falling Edges in the 
Light Curve 

The basic shape of the light curve is explained by the structure of the last closed field lines, as 
discussed in Sect. ??. In that explanation, electrons are assumed to move parallel to the field 
lines, i.e. the pitch angle <p of electrons is assumed to be 0. Since the beaming angle of the 
emitted photons is ~ 1/T ~ 10~ 7 , all the photons are considered to be emitted tangentially to 
the field lines. However, this scenario explains neither the exponential behavior of the edges nor 
the energy dependence of the rise/fall time of the edges, which are clearly visible in the observed 
data. 

Observed results can be explained by, for example, assuming that the emission from each 
field line is not beamed with an angle of ~ 1/T, but has an exponential angular distribution 
characterized by its decay constant 6 C : 

F(6) = F exp(-6/6 c ) (8.10) 

where 9 C is dependent on the photon energy. Then, the resulting light curve should be broadened 
compared with the one without a sizable emission angle , as shown in Fig. 18.91 

Here, I additionally assume that the shape of the the outer edges (the rising edge of PI and 
the falling edge of P2) reflect the angular distribution of the emission along a single field line, 
while the shape of the inner edges (the falling edge of PI and the rising edge of P2) are governed 
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Figure 8. 9: Top left: A simplified top view of a pulsar of the standard model. The emission is strongly 
beamed along the field lines. Bottom left: A simplified top view of a pulsar of the assumed model. The 
emission has an exponential angular distribution along the field lines. Top right: The emission profile map 
of the standard model assuming that inclination angle between the rotation axis and the dipole axis is 65 
degrees. Figure adopted from H53\l . Bottom right: The hypothetical light curve. A black histogram shows 
the basic structure of the light curve explained by the emission profile map shown in the top panel. Red 
lines show the hypothetical broadening effect taking into account the angular distribution of the emission 
Original figure is adopted from and red lines are added by myself. 

by the overall field line structure (see Fig. 18 .9b . From Eq. 17.41 and Eq. 17.71 the rise time of PI 
and the fall time of P2 as a function of energy are: 

r ~ (0.022 ±0.0018) - (0.0095 ± 0.0034) log 10 (£[GeV]) (8.11) 

where the simple mean of r?l e and rf^i is adopted. Since the emission angle 6 can be translated 
into the pulse phase p as p = 9/2%, Eq. 18.1 1 H eads to 

9 C ~ (0.14 ± 0.01) - (0.06 ±0.02) log 10 (£[GeV]) [rad] (8.12) 
= (7.9 ±0.6) - (3.4 ± 1.2)log 10 (£[GeV]) [deg] (8.13) 

This is by far larger than 1/T ~ 10 -7 . The large emission angle should be the consequence 
of the large pitch angle of the high energy electrons. These electrons with a spiral orbit should 
emit synchrotron radiation, which should explain the observed gamma-rays at least below 3 GeV 
(where the energy dependence of the edges is clearly seen as shown in Fig. 17.121 ) The strength 
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of the magnetic field at 1000 km from the neutron star is > 10 6 G and high energy gamma-rays 
are therefore reasonably expected (see Sect. ??). 

However, in order to realize such a large pitch angle orbit, there must exist a process which 
nearly instantly provides a large perpendicular momentum, because a gradual increase of per- 
pendicular momentum cannot occur due to the strong energy loss by the synchrotron radiation 
itself. (This is actually the reason why the pitch angle of electrons is normally considered to be 
0.) 

A possibility that the pitch angle <p can has a large value has been discussed by, for example, 
S. A. Petrova (see HI 5011 and references therein). The basic idea is as follows: The accelerated 
electrons efficiently absorb radio photons with a frequency corresponding to the cyclotron fre- 
quency when it is converted into the electron rest frame. This cyclotron resonance absorption of 
radio photons is so efficient that it can cause a large pitch angle even though electrons are contin- 
uously losing energy by the synchrotron radiation. This scenario is used in [|93l to reproduce the 
Crab pulsar energy spectrum at around 100 MeV (see Fig. I8.1i However, this process produces 
neither the exponential angular distribution nor the electron-energy dependent pitch angle (see 

C31). 
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Figure 8.10: Left: Radio delay as a function of energy. A simple linear function fit to X-ray range is 
indicated by a dashed line. Figure adopted from H138\l Right: The same as the left panel but Fermi-LAT 
and MAGIC results analyzed by myself are shown. The energy dependence calculated based on Eq. 17.751 
is shown as a solid line above 100 MeV. The dotted line is an extrapolation of it down to 1 MeV. The points 
below 1 MeV are copied from the left panel by eye. 



In Fig. 18.101 the first peak (PI) position as a function of energy is shown from optical to 
high energy gamma-rays. It is expressed as "radio delay", which tells how much earlier the 
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peak at a given energy occurs compared to the radio peak (Jodrell Bank, at 610 MHz). The 
energy dependence of the radio delay has been studied in the X-ray range. As discussed by many 
authors (see e.g. H116I . 11 1 3 VII and (11381 ). the simplest explanation for the radio delay and its 
energy dependence is the difference in the emission region within the magnetosphere. i.e., the 
higher the energy is, the more inner of the magnetosphere the emission region is. In [113811 . a 
simple linear function was fitted to the data, deriving a shift of 0.6 ± 0.2/xs keV -1 (see the left 
panel of the figure). It corresponds to 

L(E) = (180 ±60) x (E[keV}) [m] (8.14) 

where L indicate the path length difference with respect to the radio emission region. Eq. 17.151 
can be rewritten as (130 ± 20) — (70 ± 30)log 10 (£'/GeV) [fxs] in radio delay, which is shown 
as a black line in the left panel of the figure. In path length, this corresponds to 

L(E) = (40 ±6) - (20±9)log 10 (£[GeV]) [km] (8.15) 

Considering the size of the light cylinder (Rl = 1500 km), Eq. 18. 141 and Eq. I8.15l are reasonable. 
On the other hand Eq. 18.141 and Eq. 18.151 are apparently not consistent as can be seen from 
Fig. 18.101 Actually, the shape of the energy spectrum of the Crab pulsar (see Fig. ??) and 
the theoretical model (see Fig. ??) suggests different mechanisms for X-ray and high energy 
gamma-ray radiation. Therefore, it may be natural that the energy dependence of the radio delay 
is different in the two well- separated energy ranges. 

It should be noted that the definition of the peak phase is different for different types of 
analyses. For example, A. A. Abdo et al. (see flU), which is the official publication from the 
Fermi collaboration, determined the peak phase by fitting an asymmetric Lorentzian function 
to the binned light curves. In their analysis, the peak phase above 100 MeV is constant and 
—0.008 ± 0.001, which corresponds to 280 ± 30 /is in radio delay, though their highest energy 
bin of 3 to 10 GeV shows marginal deviation of —0.005 ± 0.002, corresponding to 170 ± 70 fis 
in radio delay. 
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Development of High Quantum Efficiency 
Hybrid Photodetector HPD R9792U-40 

see http://wwwmagic.mppmu.mpg.de/publications/theses/TS aito .pdf 
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Chapter 10 

Conclusions and Outlook 



10.1 Conclusions 

Before 2007, the energy spectrum of the Crab pulsar had been measured only up to ~ 5 GeV by 
the satellite-borne detector, EGRET, while IACTs had set flux upper limits only above 100 GeV. 
There existed no sufficient measurement at around the cut-off energy, i.e, at energies between a 
few GeV and a few tens of GeV, while the spectral shape at around the cut-off energy is essential 
to constraining the emission region of the pulsation, 

The MAGIC telescope with the newly implemented SUM trigger successfully detected emis- 
sion from the Crab pulsar above 25 GeV during the observations between October 2007 and 
February 2008 thanks to the collective efforts of my colleagues, T Schweizer, M. Lopez, A. N. 
Otte, M. Rissi and M. Shayduk. However, an in-depth analysis and detailed discussion in com- 
parison with the adjacent energy range had not yet been performed. Also, a new satellite-borne 
detector, Fermi-LAT, became operational in August 2008 and the observational data were made 
public in August 2009. 

In this thesis, the Crab pulsar has been studied in detail in the previously (almost) unstudied 
energy gap between 5 GeV and 100 GeV. For the analysis, I used the data from both the upgraded 
MAGIC telescope and the public Fermi-LKT data. The main results are summarized as follows: 

• MAGIC observations between October 2007 and January 2009 resulted in the detection 
of the Crab pulsar above 25 GeV with a statistical significance of 4.3 a, 7.4 a and 7.5 a 
for the first peak (PI), the second peak (P2), and the sum of the two peaks (PI + P2), 
respectively. 

• The energy spectrum of the Crab pulsar is consistent with a power law with an index 

of 3.5 ± 0.5 between 25 GeV and 100 GeV for PI, P2 and PI + P2. At 30 GeV, the 

flux of P2 is twice as high as that of PI. 

• A variation of the flux and the light curve of the Crab pulsar on a yearly time scale were 
not found in the MAGIC data. 
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• One year of Fermi-LAT data showed a clear detection (> 100 a) of the Crab pulsar from 
100 MeV to ~ 30 GeV. Between 100 MeV and ~ 30 GeV, the energy spectrum is consis- 
tent with a power law with an exponential cut-off, for total pulse (TP), PI, P2 and PI + P2. 
The cut-off energies are estimated to be 6.1 ± 0.5 GeV, 3.7 ± 0.3 GeV, 5.9 ± 0.7 GeV and 
4.5 ± 0.3 GeV for TP, PI, P2 and PI + P2, respectively. Due to the small detector area of 
~ 1 m 2 , the statistical uncertainty of the spectrum above ~ 10 GeV is rather large and it is 
not possible to detect any pulsed signal above 30 GeV. 

• From the Fermi-LAT observations, the superexponential cut-off assumption (r 2 = 2.0 
in Eq. \6M is ruled out by 4.8 a, 5.0a, 4.3a and 7.7a for TP, PI, P2 and PI + P2, 
respectively. 

• The combination of the results from MAGIC and Fermi-LAT revealed that the ex- 
ponential cut-off spectra determined by Fermi-LAT are inconsistent with MAGIC re- 
sults above 25 GeV by > 2.1a, > 4.3a and > 5.3a for PI, P2 and PI + P2, respectively, 
even if the possible absolute energy scale difference between the two experiments is 
carefully taken into account (up to 30%). 

• The flux ratio of P2 to PI and that of Bridge to PI increase rapidly with energy between 
100 MeV and 100 GeV. This behavior is similar to that in the energy range below 1 MeV 
but contrary to that in the energy range between 1 MeV to 100 MeV. 

• Both edges of the two peaks show a clear exponential behavior. In addition, the outer 
edges, i.e., the rising edge of PI and the falling edge of P2 become sharper as the 
energy increases, while the inner edges, i.e. the falling edge of PI and the rising edge 
of P2 have no energy dependence. The rise time of PI (r^ e ) and the fall time of P2 

( T faii) can De expressed as 

r£l e (E) = (2.02 ± 0.08) x 10~ 2 - (9.4 ± 1.3) x l(T 3 log 10 (£[GeV]) 
Tfau(E) = (2.42 ± 0.16) x 10~ 2 - (9.6 ± 3.1) x 10- 3 log 10 (£[GeV]) 

• The phase of the first peak has a slight but significant energy dependence. This shift 
can be expressed as 

Peakl(E) = (-3.8 ± 0.6) x 10" 3 + (2.1 ± 0.9) x 10~ 3 log 10 (£[GeV]). As the energy 
increases, the peak position shifts to a later time in the light curve. For the second peak, 
because of the broader width, the peak phase is determined with a worse precision, and no 
significant energy dependence has been found. 

• In the Fermi-LAT data above 10 GeV, a hint of the third peak is seen at phase ~ 0.75 with 
a significance of 3.5 a. However, in the MAGIC data, only a 1.7 a excess has been found 
and the flux upper limit based on the MAGIC data is in marginal contradiction with the 
Fermi-LAT results. 

• Aiming for better observations of pulsars and other sources below 100 GeV with MAGIC, 
I participated in the development of a new photodetector, the Hamamatsu hybrid pho- 
todetector HPD R9792U-40. Compared to the currently used PMTs, its photodetection 
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efficiency is twice higher and its ion-feedback rate is 500 times lower. Its charge resolu- 
tion is excellent, too. The lifetime of the photocathode was proven to be long enough to 
allow a ten year observation time without significant degradation. A compensation circuit 
for the correction of the temperature dependence of the gain and safety circuits against the 
strong light were also successfully developed. 

Based on these results, the following physics conclusions have been drawn: 

• The extension of the pulsed gamma-ray emission up to 100 GeV observed by MAGIC sets 
a lower limit in height of the emission region at 7.8 times the neutron star radius. This rules 
out the inner magnetosphere emission scenario, i.e. the Polar Cap model, for the pulsation 
mechanism. Strong rejection of the super-exponential cut-off assumption by Fermz'-LAT 
also favors the outer magnetosphere emission scenario, i.e. the Slot Gap model or the 
Outer Gap model. 

• The rejection of the exponential cut-off assumption by the combined analysis of Fermi- 
LAT results and MAGIC results requires modifications of the standard outer magneto- 
sphere model. If the magnetic field has an ideal dipole structure, there must be a place 
where the acceleration electric field is more than 10 times larger than that of the standard 
model (> 3.3 x 10 7 [V/cm]). A distorted dipole structure of the magnetic field is another 
possible explanation. 

• It is unlikely that the contribution of the inverse Compton scattering is the reason for the 
discrepancy between the standard outer magnetosphere model and the observed results, 
considering the energy of accelerated electrons and that of possible target photons. 

• The radiation efficiency above 100 MeV is estimated to be (1.36 ± 0.04) x 10 -3 from the 
Fermi-LAI measurement while that above 30 GeV is estimated to be (8.4 ± 1.1) x 10~ 6 
from the MAGIC measurement. The discrepancy in the radiation energy above 30 GeV 
between the standard outer magnetosphere model and the MAGIC measurement amounts 
to 0.6% of the radiation energy above 100 MeV. 

• The exponential behavior of the pulse edges can be explained by assuming that the emis- 
sion angle with respect to the magnetic field line has an exponential distribution. Under this 
assumption, the energy dependence of the exponential decay constant 9 C can be expressed 
as 6 C = (7.9 ± 0.6) - (3.4 ± 1.2) log 10 (£[GeV]) [deg], E < 100 GeV. 

• The simplest explanation for the energy dependence of the peak phase is that the emission 
region shifts inward toward the neutron star. The energy-dependent difference in path 
length L(E) with respect to the radio emission region can be written as 

L(E) = (40 ± 6) - (20 ± 9)log 10 (£/GeV) [km], E < 100 GeV. 
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10.2 Future Prospect: Observations of Other Gamma-ray Pul- 
sars 

The power-law-like extension of the gamma-ray energy spectrum of the Crab pulsar beyond the 
cut-off energy is a new discovery and currently there is no concrete theoretical explanation for 
it, except the discussions presented in this thesis. It would be necessary to check if this feature 
is unique for the Crab pulsar or common for all/some of the other pulsars. Although Fermi-LPR 
detected 46 gamma-ray pulsars, it is not possible to study the spectral behavior well beyond the 
cut-off energy with Fermi-LAI data due to the limited effective area of the detector. MAGIC is 
currently the only detector that can study pulsars at energies well beyond the cut-off energy. The 
next target of observation for MAGIC could be the Geminga pulsar, which is the second brightest 
above 1 GeV (next to the Crab pulsar) among the pulsars in the sky region which MAGIC can 
observe. Above 10 GeV, its flux is comparable to that of the Crab pulsar, according to the Fermi- 
LAT observations. 

10.3 Future Prospect: Improvement of the Telescope Perfor- 
mance 

The Crab pulsar is the only pulsar that has up to now been detected by an IACT Other pulsars 
have never been detected from ground despite tremendous efforts. Even for the Crab pulsar, the 
energy spectrum could be determined only with the moderate statistical significance after 59.1 
hours of observations, which is relatively long for IACT observations. The measured energy 
spectrum of the Crab pulsar is consistent with a power law above 25 GeV. However, if the sta- 
tistical uncertainty of the measurements is smaller and the energy resolution is better, a possible 
curvature of the spectrum might become visible, which surely helps to understand the reason for 
the spectral extension after the cut-off. 

In order to detect more pulsars and determine the energy spectrum with higher precision, one 
needs to meet the following requirements: 

• A lower energy threshold. 

• A larger effective area below 100 GeV 

• A better (hadron+muon)/gamma separation below 100 GeV 

• A better angular resolution below 100 GeV 

• A better energy resolution below 100 GeV 

For meeting these requirements, several improvements might be explored, as described in the 
following subsections. 
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10.3.1 Installation of the HPD R9792U-40 

The replacement of the PMTs in the camera with the HPDs described in Chapter ?? would im- 
prove the telescope performance. The HPDs will double the number of detected Cherenkov pho- 
tons from air showers. One can record shower images with higher precision, which will lead to a 
better (hadron+muon)/gamma separation, a better angular resolution and a better energy resolu- 
tion. The energy threshold will also be lowered. The effective area should also increase largely, 
not only thanks to the higher photodetection efficiency but also due to the lower ion-feedback 
rate, allowing a more efficient trigger. However, the effect of the fast-and-huge afterpulsing, 
which may be attributed to the generation of characteristic X-rays inside the HPD, should be 
carefully studied. 



10.3.2 Stereoscopic Observation with the MAGIC Stereoscopic System 




Figure 10.1: Left: A photograph of the two MAGIC telescopes adopted from H215H . Right: Concept for the 
stereo observations. Using two images recorded by the two telescopes, the arrival direction, the shower 
maximum height and the impact point of an air shower can be reconstructed with a much higher precision 
than with a single telescope. 

In October 2009, MAGIC started stereoscopic observations with two telescopes. As schemat- 
ically explained in the right panel of Fig. 110.11 if an air shower image is recorded by the two 
telescopes, the arrival direction can be reconstructed better than by a single telescope. In addi- 
tion, the impact point of the shower and the shower maximum height can be determined with 
higher precision. This leads to a better energy resolution. The information of the shower max- 
imum height has another advantage. A muon image may look like a low energy gamma-ray 
image if the impact distances from both telescopes are large (> 80 m). However, the majority of 
muons can be identified by the reconstructed shower maximum height (a muon does not produce 
a shower but can create a shower-like image). In the case of muon images, the height should 
be reconstructed to be ~ 5 km, which is unusually low for low energy gamma-rays (see Fig. 
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110.21) . Also, Hadron/gamma-ray separation will improve simply due to the double amount of 
information from the two telescopes. 



Shower Maximum Height: Size < 100 | | Shower Maximum Height: 100 < Size < 500 | | Shower Maximum Height: Size > 500 




Shower Maximum Height [km] Shower Maximum Height [km] Shower Maximum Height [km] 



Figure 10.2: The shower maximum height distributions for observed data (black histograms) and gamma- 
ray MC (red histograms). SIZE < 100 (left), 100 < SIZE < 500 (middle) and SIZE > 500 (right) 
are shown. At SIZE > 500, two distributions are not clearly separated. At 100 < SIZE < 500, two 
peaks can be seen for the observed data. The first peak at ~ 5 km can be explained by the the contribution 
from large-impact-distance muons (see text for detailed explanation). At SIZE < 100, only one peak 
can be seen at ~ 6 km for the observed data, suggesting that most of the events are large-impact-distance 
muons. Consequently, the distribution of the observed data is clearly separated from that of the gamma- 
ray MC. Therefore, the shower maximum height, which can precisely be estimated only with the stereo 
observations, is a powerfid parameter for the gamma- ray/liadron separation for very low energies. 

However, the requirement for a coinciding signal from the two telescopes reduces the effec- 
tive area especially for low energies. Below 50 GeV, the reduction can be as large as a factor of 
~ 5. The energy threshold also increases. The trigger condition, i.e. whether or not the coinci- 
dence condition is fulfilled, must be carefully studied taking into account the advantage and the 
disadvantage of stereoscopic observation. 



10.3.3 Pulsar Observation with CTA 

A new project in VHE gamma-ray astronomy named Cherenkov Telescope Array (CTA) was 
recently proposed and a large international collaboration was formed (see H20910 . The basic 
concept of CTA is an array of a large number of IACTs spread over a square kilometer area. 
The telescope configuration would be as follows; 

• The array consists of 3 different sizes of telescopes, 3 — 5 large size ones, 20 — 30 middle 
size ones, and 20 — 30 small size ones. 

• The large size telescopes have a reflector diameter of ~ 23 m with a FoV of 5 degree. 

• The middle size telescopes have a reflector diameter of ~ 12 m with a FoV of 8 degree. 
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Figure 10.3: Left:An illustration of the CTA project. Figure adopted from H209\l . Right: Sensitivity curves 
for currently operational experiments (GLAST(Fermi-LAT), MAGIC, HESS, VERITAS and Milagro) and 
expected sensitivity curves for future experiments (HAWC and CTA). Figure adopted from H2II\I . 

• The small size telescopes have a reflector diameter of ~ 6.5 m with a FoV of 10 degree. 

• The distance between telescopes is around the order of 100 m. 

• Two arrays will be constructed, a larger one in the southern hemisphere and a smaller one 
comprising only 23 m and 12 m telescopes in the northern hemisphere. 

The parameters described above such as the numbers of telescopes, the diameters of reflectors 
and the angle of FoVs are not yet fixed. The primary purpose of the CTA project is the obser- 
vation of VHE gamma-ray sources in the energy range between 100 GeV and 10 TeV with a 
sensitivity 5 — 10 times higher than that of current IACTs. The accessible energy range will also 
be extended down to 10 - 20 GeV and up to 100 TeV. 

For pulsar observations, a sensitivity below 100 GeV is essential. Below 100 GeV, the sensi- 
tivity is mainly determined by the large size telescopes, whose reflector area will be ~ 1.7 times 
larger than that of MAGIC. If the SUM trigger system is used in CTA, and if a coincidence be- 
tween multiple telescopes is not required for a trigger, then the energy threshold of the CTA may 
be roughly estimated to be ~ 15 GeV by scaling with the reflector diameter. Since there will be 
3-5 telescopes, the effective area would increase correspondingly. The coincidence requirement 
might increase the threshold energy, but the reduction of the effective area may not be as much 
as in the case of the MAGIC-stereo system. The gain in sensitivity below 100 GeV should be 
intensely studied before all the designs are fixed and the construction of telescopes starts. 
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